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Abstract: 

In  this  paper  we  present  a  methodology  for  checking  the  local  quality  of  recipes 
for  the  recovery  of  stresses  or  derivatives  from  finite  element  solutions  of  linear 
elliptic  problems.  The  methodology  accounts  precisely  for  the  factors  which  affect 
the  local  quality  of  the  recovered  quantities,  namely,  the  geometry  of  the  grid,  the 
polynomial  degree  and  the  type  of  the  elements,  the  coefficients  of  the  differen¬ 
tial  equation  and  the  class  of  solutions  of  interest.  We  give  examples  of  how  the 
methodology  can  be  used  to  obtain  precise  conclusions  about  the  quality  of  a  class 
of  recipes,  based  on  least-squares  patch-recovery,  in  the  interior  of  complex  grids, 
like  the  ones  employed  in  engineering  computations.  By  using  this  approach  we 
were  able  to  discover  recipes  which  are  much  more  robust  than  the  ones  which  are 
currently  in  use  in  the  various  finite  element  codes. 
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1  Introduction 


The  subject  of  superconvergence  studies  techniques  which  can  be  employed  to 
increase,  at  a  relatively  small  cost,  the  accuracy  of  quantities  (derivatives  of  the 
solution,  strains,  stresses  etc.)  recovered  from  finite  element  solutions.  Because 
of  its  obvious  importance  in  finite  element  practice  the  superconvergence  of  finite 
element  approximations  has  been  studied  by  several  investigators  in  the  mathe¬ 
matical  and  the  engineering  literature;  see  for  example  [1-46]  and  the  citations  in 
these  papers. 

We  distinguish  two  types  of  superconvergence  according  to  the  way  that  the 
values  of  the  solution  quantities  are  recovered,  namely: 

1.  Direct  superconvergence:  The  values  of  the  solution  quantities  are  obtained 
by  direct  evaluation  from  the  finite  element  solution  at  special  points ,  the 
superconvergence  points  for  the  solution  quantity  in  the  element  (see,  among 
others,  [6],  [10],  [13],  [14],  [47]  and  [50]). 

2.  Superconvergence  via  averaging:  The  values  of  the  solution  quantities  at  a 
point  are  obtained  by  averaging  the  values  of  the  finite  element  solution  in  a 
neighborhood  of  the  point.  According  to  the  definition  of  the  neighborhood, 
in  which  the  averaging  is  computed,  we  have: 

a.  Global  or  subdomain  averaging:  The  neighborhood  is  either  the  entire 
domain  (see  [2],  [3],  [9],  [11],  [12],  [15],  [16-17],  [42])  or  a  subdomain 
which  does  not  depend  on  the  mesh-size  (see  [15],  [16-17],  [42]). 

b.  Local  or  semilocal  averaging:  The  neighborhood  is  a  mesh-cell  of  ele¬ 
ments  which  are  connected  to  a  vertex,  an  edge  or  an  element  (local 
^-neighborhood)  or  such  a  mesh-cell  and  a  few  mesh-layers  around  it 
(semilocal  ^-neighborhood)  (see  also  [46]). 

In  this  paper  we  will  study  superconvergence  via  local  averaging;  the  case  of 
direct  superconvergence  will  also  be  considered  (as  a  special  case  where  the  aver¬ 
aging  is  the  pointwise  evaluation  of  the  solution  quantity  directly  from  the  finite 
element  solution).  For  a  detailed  study  of  direct  superconvergence  for  meshes  of 
triangles  and  quadrilaterals  see  [47-50]. 

The  objective  of  this  paper  is  to  develop  computer-based  methodologies  for 
checking  the  quality  of  recipes  for  the  recovery  of  solution  quantifies.  We  will  apply 
the  methodologies  in  the  study  of  recipes  based  on  least-squares  patch-recovery 
of  the  type  given  in  [35-40];  these  recipes  recover  the  solution  quantities  in  an 
element  from  the  finite  element  solution  in  the  element  and  one  mesh-layer  around 
it.  It  is  also  possible  to  define  various  recipes  (see  [41]  and  [44]  and  Section  4 
below)  which  take  into  account  the  information  about  the  data  (the  body-force, 
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boundary-traction)  and  the  differential  operators  (i.e.  the  differential  equation,  the 
traction  boundary-condition)  and  the  quality  of  these  recipes  can  be  checked  by 
using  exactly  the  same  approach. 

In  the  engineering  literature  the  quality  of  recipes  is  often  evaluated  by  observ¬ 
ing  their  performance  in  a  few  example  problems  ( benchmarks )  (e.g.  [24],  [35-40], 
[41],  [43]  and  [44]).  This  approach  has  serious  shortcomings  because  the  quality 
of  recipes  is  sensitive  to  the  structure  of  the  solution  and  the  meshes  employed. 
The  quality  of  the  recovered  solution  quantities  can  be  very  different  depending  on 
whether  the  meshes  are  translation-invariant  or  general,  whether  the  meshes  are 
nearly  equilibrated  (or  not)  and  whether  the  solution  has  singular  points.  The  per¬ 
formance  of  a  recipe  can  be  very  different  in  the  neighborhood  of  a  singular-point, 
at  a  smooth- boundary  and  in  the  interior  of  the  mesh.  Hence  benchmark  computa¬ 
tions  could  motivate  misleading  conclusions  unless  each  benchmark  is  more  or  less 
representative  of  a  precise  class  of  computational  problems  and  the  conclusions  are 
employed  for  this  class  of  problems. 

In  this  paper  we  will  define  precisely  the  notion  of  the  quality  of  the  recipes 
and  we  will  present  a  theory-based  computational  methodology  for  the  assesment 
of  the  quality  of  the  recipes  in  the  interior  of  the  domain  and  smooth  solutions. 
Similar  approaches  for  checking  the  quality  of  recipes  at  smooth  boundaries  and 
in  the  neighborhood  of  singular  points  will  be  given  in  forthcoming  papers.  The 
approach  takes  directly  into  account  the  factors  which  affect  the  quality  of  the 
recipes,  namely: 

a.  The  local  geometry  of  the  grid ; 

b.  The  coefficients  of  the  differential  equation ; 

c.  The  polynomial  degree  and  the  type  of  the  elements  (triangles,  tensor-product 
quadrilaterals,  serendipity  squares  etc.); 

d.  The  class  of  rotations  of  interest. 

Below  we  will  give  illustrative  examples  of  how  the  quality  of  a  recipe  is  influenced 
by  these  factors. 

Following  this  Introduction  we  outline  the  model  elliptic  problems  (Poisson’s 
equation  and  the  elasticity  problem).  We  give  a  new  definition  of  superconvergence 
and  based  on  this  definition  we  develop  a  theory  for  checking  recipes  in  the  interior 
of  any  class  of  meshes  and  for  any  class  of  material  coefficients  and  we  give  several 
examples  of  application  of  the  theory  for  studying  the  quality  of  recipes  based  on 
least-squares  patch-recovery  in  the  interior  of  complex  finite  element  grids.  We 
give  examples  of  new  recipes  (which  are  similar,  but  less  expensive,  to  the  ones 
given  in  [35-40])  which  also  employ  the  available  information  about  the  differential 
equation  and  the  right-hand  side  in  the  design  of  the  recipe.  By  employing  the 
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methodology  of  the  paper  we  were  able  to  show  that  one  of  these  recipes  is  much 
more  robust  than  the  recipes  in  [35-40]. 

2  The  model  problems 

We  shall  consider  the  vector-valued  boundary- value  problem 

£.(u)  :=  -  53  lT~{aii(u))  =  fi  in  Q 

j= i  oxi 

u,  =  0  on  Tp  (1) 

2 

£V,j(u)n_,  =  U  on  Tat 

i=i 

where  *  =  1,  2. 

Here  0  C  is  a  bounded  domain  with  boundary  dQ  =  To  U  IV  ; 

n  (ni,n2)  is  the  outward  pointing  unit-normal  on  IV; 

fi,  i  —  1,2  are  the  components  of  the  load-vector  ( body-force ); 

ii ,  t  s  1,2  are  the  components  of  the  normal-flux  vector  ( traction )  applied  on  IV; 

IV  7^  0,  IV  H  IV  =  0;  u  =  (tii,u2)  is  the  solution-vector  ( displacement ); 

Ci,(u)  :=  Kf*”  +  »  *«(») :=  £  .  *.J  =  1,2  (2) 

are  the  components  of  the  flux  (strain,  stress ); 

aijkt ,  i,j,k,£  =1,2,  are  the  material-coefficients  ( elastic  constants)  which  in  the 
case  of  isotropic  plane  elasticity  are  given  by  a,j«  =  fi(8tJ8kt  +  8^8 kj )  +  A 
where  ^  is  Kronecker’s  delta  and  A ,  p  are  Lame’s  constants. 

We  also  introduce  the  scalar  elliptic  boundary- value  problem  ( heat-conduction 
in  orthotropic  medium),  namely 
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1,w 

in  fl 

u  =  0 

on  Td 

2 

Y  9*(u)n*  =  9 

k=l 

on 

(3) 


Here  /  is  the  source  function  ( heat-source );  g  the  boundary-flux  ( heat-flux );  qk(u)  := 
2  dtA 

£  Kkt  q — »  Jfc  =  1,2  are  the  components  of  the  flux- vector  ( heat-flux );  Kkt,  k,t  = 

1,2,  are  the  entries  of  the  thermal-conductivity  matrix  which  is  symmetric,  positive 
definite.  Below  we  will  let  denote  the  principal  values  of  the  thermal- 

conductivity  matrix. 

Let  us  now  cast  the  problems  in  the  variational  form.  Let  us  denote  the  space 
of  test-functions  by 


rD  :=  { («i, ”2)  |  ^  €  Hl( fl) ,  Vi  =  0  on  rD| 


H 

The  variational  form  of  the  boundary- value  problem  (1)  is  now  posed  as: 
Find  u  €  Hpo  such  that 

Bn(u,  v)  =  /  £  fiVi  +  I  Y  v  v  ^  HrD 
Ja  ,=1  J rN  ,=1 

where  the  bilinear  form  Bn  :  Hpfl  x  HpD  -»  R  is  defined  by 

dttj  dv{ 


Bn(u,v):=J  £ 


2 

*  wu:  l/l 

*  a'jkt8x  dl 

The  energy-norm  over  any  subdomain  5  C  fl  is  defined  by 

HMDs  :=  y/Bs(v,v) 


(4) 


(5) 


(6) 


(7) 


where  Bs(u,v)  has  the  obvious  meaning. 

In  the  case  of  the  scalar  elliptic  problem  (3)  the  bilinear  form  is  given  by 

6o(ti,t;)  :=  f  Y  Kktir~  "3~-  The  weak-solution  of  (3)  satisfies: 

Jn k,=1  aXf  oxk 
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(8) 


Find  u  €  HlD  :=  |v  €  Hl(Q)  |  t;  =  0  on  To  |  such  that 

M«>v)  =  Lfv  +  L  iv  V  Hi 

Jo  Jr  s 

The  energy-norm  in  any  subdomain  S  C  (l  is  defined  by  |||u|||s  :=  ^bs(v,v)  . 

Let  T  =  {T/,}  be  a  family  of  meshes  of  triangles  or  quadrilaterals  with  straight 
edges.  It  is  assumed  that  the  family  is  regular  (for  the  triangles  the  minimal  angle 
of  all  the  triangles  is  bounded  below  by  a  positive  constant,  the  same  for  all  the 
meshes;  for  the  quadrilaterals  see  conditions  (37.40)  in  Ciarlet  [52]).  The  meshes 
are  not  assumed  to  be  quasiuniform.  We  introduce  the  conforming  finite-element 
spaces  for  the  scalar  model  problem  (the  corresponding  spaces  of  vector-valued 
functions  for  the  elasticity  problem  are  defined  similarly): 

SJ(T,)  :=  {u  €  C°(n)  |  «|T,  €  Sf(r4) ,  k  =  1,. . .  ,Af(T,)}  ,  (9) 

where  5J(rt)  =  VArk),  for  the  meshes  of  triangles,  while  for  the  meshes  of  quadri¬ 
laterals 


SIM  ■■=  {-»  €  C~(r,)|  »oF„e  (10) 

where  f  :=  ( — 1,  l)2  is  the  master-element  and  FTfc  is  the  bilinear  mapping  of  f 
onto  rk;  M(Th)  is  the  number  of  elements  in  Th  and 

[f>\  P(i„i,)=  £  Oyijii}  (11) 

0<tj<p 


Let  P(i„i2)  =  £  W  x™x2  be  a  polynomial  of  degree  p;  noting  that 

0<m+n<p 


X1  =  S  °ij  X’l  X2  » 

ij=0 


=  E  h„  i\  H 

>J=0 


(12) 


we  get 


P(*i,x2) 


0<m+n<p  ij= 0  k,l= 0 


j  **  - 1 

XX  x2 


0<*<p 

0<«<p 


(13) 


Hence  any  polynomial  of  degree  p  belongs  to  the  span  of  the  shape-functions  5^(Tfc ) 
on  every  physical  quadrilateral  element  rfc. 
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We  let  S*Xd  :=  5l(rk)n^a  denote  the  discrete  test-space.  The  finite  element 
solution  uK  (for  the  heat-conduction  problem)  satisfies: 

Find  uk  G  5£rD  suc^ 

=  Lfvk+  L  ^vh  v  ”h  €  r  (14) 

^n  JrN  u 

The  error  is  eA  :=  u  —  uA.  The  finite  element  solution  and  the  error  for  the  elasticity 
problem  are  defined  similarly. 

3  A  new  definition  of  superconvergence 

The  assessment  of  the  quality  of  the  recipes  will  be  based  on  the  following 
new  definition  of  superconvergence  (see  also  [47-50]):  Let  {uh}  be  a  one  parameter 
sequence  of  finite  element  solutions  of  a  problem  which  are  computed  using  a 
sequence  of  meshes  T  =  {Th}  and  let  u  denote  the  exact  solution.  Let  us  assume 
that  we  are  interested  in  the  values  of  the  linear  functional  F(u)(x)  (i.e.  the  values 
of  the  components  of  strain,  stress  or  the  gradient  of  the  displacement)  at  a  point 
x.  We  will  let  F(ufc)(x)  to  denote  the  recovered  value  of  F(u)(x)  obtained  from 
the  averaging  T  which  is  assumed  to  be  in  the  form 

^K)(x)  :=  I>f(u)(xK(Xi)  (15) 

«=i 

where  {x^  is  a  set  of  points,  which  may  include  x,  in  an  ^-neighborhood  of 
x  (as  shown  in  Fig.  1)  and  are  the  coefficients  corresponding  to  the 

quantity  F(u).  If  F(u)  has  several  components  it  is  understood  that  (15)  holds 
with  a  different  set  of  coefficients  for  each  component.  Let  us  denote  by 

:=  (J  t  ,  r  €  Th  (16) 

S€f 

the  mesh-cell  which  consists  of  the  elements  which  include  x  in  their  interior  or 
boundary.  (If  x  is  at  a  vertex  (resp.  on  an  edge)  u;s  is  the  set  of  elements  connected 
to  the  vertex  (resp.  to  the  edge)).  We  let 

tf(u>g)  :=  max  |  F(u  -  ufc)(x)  |  (17) 

Here  J  •  |  indicates  the  modulus  employed  for  the  quantity  F(u).  (If  F(u)  is  a 
scalar-quantity  (one  component  of  stress,  strain  or  of  the  gradient)  |  •  |  is  the 
absolute  value,  if  F(u)  is  a  vector- valued  quantity  J  •  J  is  the  vector-modulus, 
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if  F(u)  is  a  second-order  tensor  J  •  |  denotes  an  invariant,  like  for  example  the 
second  invariant,  the  second  invariant  of  the  deviator,  etc.) 

We  will  be  interested  in  the  values  of  the  error  in  the  recovered  quantity  ^‘(ufc)(x) 
relative  to  the  maximum  error  of  the  finite- element  solution  in  o/s,  namely 


[  F(u)(x)-: 

e(x;F;f;u;u,„Th):= 


-  ^M(x) 


) 


0 


,  if  ¥(ws)  ±  0 


,  if  =  0 


(18) 


Let  us  note  that  the  relative  error  with  respect  to  ^(w*)  has  clear  meaning.  We 
compare  the  quality  of  the  recovered  value  at  x  with  the  maximum  error  in  the 
quantity  computed  directly  from  the  finite  element  solution  in  the  mesh-cell  u>s. 

If  the  point  x  and  the  averaging  T  are  such  that 


Q(x;F\f;u;uiyTh)  <  ^  as  h  — ►  0  (19) 

then  F(uh)(x)  will  be  called  u-ri%-superconvergent  relative  to  the  exact  solution  u 
and  the  family  of  meshes  T.  Given  a  class  of  solutions  U,  we  let 


'H{x\F\F\U\ui,Th)-  mgxQ(x,F\?,u-, ut,Th)  (20) 

denote  the  maximum  relative  error  in  F(uh)  at  x  for  the  class  U.  The  dependence 
of  H(x",F\f;ll-,ijJx,Th)  on  U  is  very  essential.  In  practical  computation  we  know 
a-priori  some  of  the  properties  of  the  solution.  For  example  if  the  distributed  load 
is  zero  we  know  that  the  solution  satisfies  the  homogeneous  differential  equation. 
Then  we  can  restrict  the  set  U  to  the  set  of  solutions  of  the  homogeneous  differential 
equation.  (Below  we  will  see  that  the  quality  of  a  recipe  can  be  very  different  for 
different  classes  U.  Of  course  it  is  preferable  that  the  quality  of  a  recipe  is  good 
for  a  large  class  U\  then  we  can  say  that  the  recipe  is  robust.) 

We  will  say  that  the  value  ^(ufc)(x)  is  U-rj%-superconvergent  if 

.  'H(x;F;F,U;wSlTh)  <  ^  as  h  — ►  0  (21) 

If,  in  addition,  we  have 


'H{x-,F;f-,U-,u>i,Th)  — »  0  as  h  — >  0  (22) 

then  ^(ufc)(x)  is  W-superconvergent  in  the  classical  sense.  It  means  that  the  rate 
of  convergence  of  ^(uh)(x)  to  the  exact  value  F(u)(x )  is  higher  than  the  global 
rate  of  convergence  of  the  quantity  F(uh). 
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It  is  well  known  that  classical  superconvergence  is  very  sensitive  to  the  geometry 
of  the  grid,  the  class  of  solutions  of  interest  and  the  coefficients  of  the  differential 
operator  (see  for  example  [48]  and  [50]  for  direct  superconvergence;  the  situation 
is  very  similar  for  superconvergence  via  averaging).  Let  us  note  that  if  7  in  (21) 
is  very  small  (smaller  than  a  tolerance)  we  can  still  practically  speak  about  a 
superconvergent  value  (with  a  given  tolerance),  which  is  the  reason  for  introducing 
the  notion  of  7%-superconvergence.  Note  that  7%-superconvergence,  as  classical 
superconvergence,  has  asymptotic  character  for  h  — ♦  0.  However,  unlike  classical 
superconvergence,  which  holds  only  for  special  grids  and  a  fixed  set  of  coefficients 
of  the  differential  equation,  tj  %-superconvergence  is  well-defined  for  any  classes  of 
grids  and  for  an  entire  class  of  coefficients. 

Remark  3.1.  Note  also  that  the  notion  of  superconvergence  can  be  understood  in 
different  ways.  For  example  we  can  be  interested  in  the  quantity. 


{r},rj 


e  e  |f(u)(x,)-^K)(*r)ry 

>€{r}ST€{»r}  ' 


(E  (*(’■))•) 1 

r€{r) 


(22) 


where  {r}  is  a  set  of  elements,  {xT}  is  a  set  of  points  in  each  element.  Then  the 
pointwise  superconvergence  is  related  to  q  =  00.  Very  often  (see  [49]  for  concrete 
examples  in  the  case  of  meshes  with  local  refinements)  we  get  — ►  0  as  h  — ►  0 

(in  fact  ~  h«).  Various  authors  are  also  calling  this  effect  superconvergence 
(see  e.g.  [13,  14]).  The  notion  of  superconvergence  depends  on  its  measure,  e.g.  if 
measured  in  Z’-norm,  q  <  00,  we  may  get  superconvergence  effects  where  for 
q  =  00  there  may  be  no  superconvergence  (see  [49]  for  further  details).  Hence  the 
term  superconvergence  has  to  be  properly  defined. 

Given  a  recipe  T  and  0  <  7  <  00  we  will  be  interested  in  the  7 %-sy.perconvergence 
regions  of  !F(uh)  in  the  element  t  £  Th  for  the  class  of  exact  solutions  U: 

.-=  (x6r|  H(x;F-,r-M-,T,Th)  <  (23) 

Given  two  recipes  Tx  and  F2  and  a  point  x  we  define  the  index  of  relative  robustness 
of  the  recipes  at  x. 


vF(u){y.-,u-,r\r2,r,Th) 


H{x\F;^U;rtTh) 

H(x;F;^;W;r,T,) 


(24) 


Here  we  will  assume  that  the  denominator  can  be  equal  to  zero  only  at  isolated 
points  (this  assumption  that  these  points  are  isolated  may  be  violated  only  for 
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very  special  mesh-geometries  like  meshes  of  squares,  triangles  arranged  in  a  three- 
directional  pattern,  etc.)  where  the  quotient  is  defined  equal  to  oo  (unless  it  is 

equal  to  —  in  which  case  we  will  define  it  by  continuity).  We  will  be  interested  in 
the  regions  in  the  element  r  €  Th  where  T2  is  more  robust  than  Tx ,  namely 

AF(umr\r'-tT,Th)  :=  {x  €  T  |  PF(u)(x;W;^l,^2;r,rfc)  >  l}  (25) 

Obviously  the  above  introduced  notions  are  defining  exactly  the  (asymptotic) 
quality  of  the  recipe.  Comparisons  of  recipes  should  be  based  on  computations  of 
the  above  quantities,  for  the  classes  of  meshes  and  solutions  of  interest,  and  not  on 
arbitrarily  selected  benchmarks. 

One  is  interested  to  know  a-priori  the  quantities  defined  above  for  classes  of 
solutions  of  interest.  (In  the  practical  computations  in  plane  elasticity  and  heat- 
conduction  the  class  of  solutions  of  interest  is  the  class  <  /  “harmonic ”  solutions 
(i.e.  the  solutions  which  satisfy  the  homogeneous  differential  equation)  with  a 
finite  number  of  algebraic  point  singularities  of  the  type  ra  at  isolated  points  on 
the  boundary  and  material-interfaces).  In  this  paper  we  will  determine  a-priori  the 
quantities  defined  in  (18),  (20),  (23),  (24),  (25)  under  the  following  assumptions: 

1.  The  grid  is  locally  periodic  with  period  h  (details  about  how  the  limit  h  — ►  0 
is  taken  are  given  below); 

2.  We  will  study  the  superconvergence  in  mesh-patches  in  the  interior  of  the 
domain  where  the  solution  is  smooth  (It  is  approximated  well  by  its  local 
Taylor-series  expansion  of  degree  (p+  1)  where  p  is  the  polynomial  degree  of 
the  elements); 

3.  The  global  modes  of  the  error  (pollutions),  in  the  elements  of  interest,  are 
negligible  when  compared  with  the  magnitude  of  the  error  in  the  local  best- 
approximation  (see  [51]  and  [53]). 

The  detailed  theoretical  setting  and  the  precise  assumptions  are  given  in  [47] 
(see  also  the  outline  given  below).  Our  theory  has  an  asymptotic  character  i.e.  it 
gives  the  precise  results  when  the  mesh  is  sufficiently  refined  or  the  accuracy  of 
the  finite  element  solution  is  reasonably  high.  The  theory  of  this  paper  cannot  be 
applied  when  the  exact  solution  is  unsmooth  (in  the  mesh-cell  or  element  of  interest) 
or  when  the  mesh  is  coarse. 

We  now  give  an  example  to  explain  the  meaning  of  the  new  definition  of  su¬ 
perconvergence.  Let  us  consider  the  Dirichlet  boundary-value  problem  for  the 
Laplacian 
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Au  =0 


in  fi  :=  (0, 1) 


u  =  u  on  3fl  j 


where  ti  is  the  boundary- value  corresponding  to  the  exact  solution  tt^,  x2)  = 
8in(TZj)  sinh(irx2).  We  computed  a  sequence  of  finite  element  solutions  {ufc}  using 
a  sequence  of  meshes  of  skewed  bilinear  quadrilaterals  }.  These  meshes  were 
obtained  by  subdividing  the  square-domain  fi  into  a  uniform  three-directional  mesh 
of  triangles  and  by  converting  each  triangle  into  three  quadrilaterals  by  connecting 
the  centroid  of  the  triangle  with  the  midpoints  of  its  edges.  The  reason  for  choosing 
such  a  mesh  instead  of  a  regular  mesh  of  bilinear  squares  is  because  we  would  like 
to  give  an  illustrative  example  of  the  new  definition  in  a  case  where  it  is  different 
than  the  classical  definition  of  superconvergence. 

We  considered  the  sequence  of  meshes  {Th},  h  =  ;  the  mesh  7\ 

,  4  o  lu  • 

is  shown  in  Fig.  2a.  We  considered  the  mesh-cell  u2h  =  (0.25  —  h,  0.75  —  h)x 

(0.25  +  h,  0.75  +  h)  which  consists  of  four  square  subcells  uj'h,  i  =  1, ...  ,4,  as  shown 

in  Fig.  2b.  We  identified  points  x**,  K  =  1, ...  ,6,  in  each  subcell  u>\.  The  points 
with  identical  subscripts  K  have  identical  relative  positions  with  respect  to  the 
corresponding  subcell.  We  computed  the  value  of 


K  =  1,.  ,,6 


where  T  is  the  ZZ-1  recipe  given  in  Section  8,  t'k  is  the  element  where  the  point 

belongs.  In  Fig.  3  we  show  the  graphs  of  0^  versus  -r  for  K  —  1,  2,  5,  6.  These 

ft 

values  converge,  for  h  — ►  0,  to  the  asymptotic  values  which  were  determined  by 
the  theoretical  analysis  of  Section  5. 

We  also  computed  the  rate  of  convergence  for  all  the  points.  From  theorem  2 
and  Remark  5.2  in  Section  5  below  we  have 


|  F(«)~  F(uk)  I  =eK*(U4)  +  c1A>+"i+«(/.'Mi),  K>0 

1=1 


where  0*  denotes  the  asymptotic  value  of  the  relative  error  at  the  point  which 

4 

corresponds  to  K  and  ’P ((J  u?J,)  is  the  maximum  error  in  F(u)  in  the  assembly  of 

i=i 

4 

the  four  mesh-cells.  We  have  'JKlJ  wJJ  =  Coh  +  o(/i);  hence  we  expect  the  error 

i=i 
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in  the  recovered  quantity  to  converge  to  zero  linearly,  unless  is  zero  (i.e.  the 
point  is  a  0%-superconvergence  point). 

In  Table  1  we  give  the  values  of  the  error  in  the  recovered  quantity  max 

i=l,.  ...4 

|  (-F(u)  —  (x*)  |  versus  h.  For  the  points  with  K  =  1,  2,  3,  5,  6  we 

employed  T  =  Fzz~x  while  for  the  point  with  K  =  4  we  let  T  =  J~DlR.  We  note 
that  we  have  ©3  =  0%  and  the  values  of  the  error  in  the  recovered  quantity  at  the 
comers  of  the  subcells  Me  superconvergent  with  the  rate  two  (i.e.  <r'Q  =  1).  On  the 
other  hand  at  the  points  which  correspond  to  K  =  1,  2  the  error  in  the  recovered 
quantity  converges  to  zero  with  the  rate  one.  For  K  =  4  it  appears  that  the  rate 
is  two.  This  is  because  0^  =  0.0162%  and  QKC0  =  7.128E-03  while  Cx  —  4.750 
and  hence  the  true  asymptotic  rate  of  the  error  in  the  recovered  quantity  cannot 
be  seen  for  the  mesh-size  considered  here.  We  have 

max,  |  (F(u)  -  *K))<*?)  |  =  +  4.75)  h? 

We  can  detect  the  correct  order  of  convergence  when  h  is  sufficiently  small  so  that 
the  first  term  in  the  parenthesis  becomes  of  the  same  order  or  overcomes  the  second 

term;  this  is  possible  for  h  <  ^  Fig.  S  shows  very  clearly  the  preasymptotic 

phase  (for  h  large),  where  our  theory  is  not  applicable,  and  the  asymptotic  phase 
(for  h  sufficiently  small),  where  our  theory  predicts  very  well  the  pointwise  values 
of  the  relative  error. 

4  Recipes  based  on  least-squares  patch-recovery 

Below  we  describe  the  ZZ-recipes  which  are  now  commonly  used  in  finite  element 
codes.  We  also  describe  a  new  class  of  recipes  which  are  similar  to  the  ZZ-recipes 
and  also  take  into  account  the  information  about  the  differential  operator  and  the 
source  term. 


4.1  Definition  of  the  ZZ-recipes 

We  will  now  outline  the  definition  of  the  ZZ-recipes  for  the  elasticity  problem. 
The  ZZ-recipes  for  the  orthotropic  heat-conduction  problem  are  defined  analo¬ 
gously  (see  also  [36-40]  for  further  details). 

Let  u>x  U  T>  denote  the  patch  of  elements  connected  to  vertex  X,  as 

X€N(  t') 

shown  for  example  in  Figs.  4a,  4b.  For  each  patch  u>x  we  recover  the  patch- 
projection  o'zzi  by  solving  the  following  problems: 

a.  Discrete  least-squares : 


13 


IkK)  "  ^2zlb(«x).c->,{y<}^  =  |k(ufc)-<rA||L^x).c->,{y^  (26) 


*5€Pp(«X) 
•J-1,2 


where  {y<}^  denotes  a  set  of  sampling  points  in  w*  and 


nap 

lklli,a(a>x).c-i,{yt}?*; :=  Y 

m=l 

b.  Continuous  least-squares: 


Y  ^i(yTO)  *«( ym) 


(27) 


lk(u»)  -  lk(“»)  -  **  llin»x).=-  <28) 

where 

ll^lll*(u,x),c-»  :=  /  Y  °ijciikt°kt  (29) 

Jw*  ij,k,t=i 

In  [36]  cijkt  was  chosen  to  be  the  material-tensor  aijkt.  It  is  also  possible  to  define 
the  patch-projections  by  letting  ciikt  =  Sik6jt. 

A  continuous  piecewise  polynomial  recovered  stress  trzz  is  obtained  over  each 
element  by  combining  the  patch-projections  axz  which  correspond  to  the  vertices 
of  the  element.  In  the  examples  the  quantity  trzz  was  constructed  as  follows: 

a.  Linear  (or  bilinear)  elements  (p  =  1). 

The  recovered  C°-continuous  piecewise-linear  (or  bilinear)  flux-field  <rzz  over 
the  domain  is  constructed  as 

<rZZ(x)  =  Y**\X  *x  (x)  (30) 

Here  crx  is  the  value  of  the  crx  at  the  vertex  X  and  is  the  piecewise-linear 
(or  bilinear)  basis  function  associated  with  this  vertex.  Here  and  below  we  omitted 
the  subscript  ZZ  from  crzz. 
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b.  Quadratic  (or  biquadratic)  elements  (p  =  2). 

For  the  quadratic  triangles  a  piecewise  quadratic  Langrangian  representation 
of  <rzz  is  constructed  from 


*ZZ(X)  =  £®-X|X  +  +  *Xi\Y)  Mx)  (31) 

where  Y  denotes  the  midside  node  for  the  edge  X1X2;  $Xi  $y  are  the  Langrangian 
quadratic  basis  functions  associated  with  the  nodes  X ,  V,  respectively.  For  the 
biquadratic  quadrilaterals  we  let 


where  X,  Y,  $y  have  the  same  meaning  as  before  and  Z  and  $z  denote  the  in¬ 
ternal  node  and  the  corresponding  shape-function,  respectively,  for  the  biquadratic 
element. 

The  quantity  trzz  for  cubic  elements  is  defined  similarly. 

Remark  4-1.  For  the  linear  triangles  we  employed  sampling  points  located  at  the 
centroids  of  the  elements  while  for  the  quadratic  triangles  the  sampling  points  are 
taken  at  the  midpoints  of  the  sides.  These  points  are  the  same  as  the  ones  given 
in  [36-39].  For  the  quadrilaterals  we  employed  the  mapped  p  x  p  Gauss- Legendre 
points  as  sampling  points  for  the  discrete  recipe;  for  the  continuous  version  of  the 
recipe  we  employed  the  mapped  10  x  10  Gauss-Legendre  quadrature  rule. 


4.2  Examples  of  “statically  admissible”  recipes 

Here  we  give  examples  of  alternative  recipes  which,  like  the  ZZ-recipes,  are 
based  on  local  least-squares  and  also  take  into  account  the  knowledge  about  the 
differential  operator  and  the  source  term.  Let  us  assume  that  we  are  recovering  the 
stresses  from  a  finite  element  solution  of  the  elasticity  problem  with  a  given  body 
force  f.  For  each  patch  u>x  we  recover  a  “ statically  admissible ”  patch-projection 
&x  by  solving  one  of  the  following  problems: 

a.  Least-squares  fit  for  the  stresses : 

Find  €  SSA  such  that 


||<r(ufc)  -  d-fj|i2(Wx)c-My<>n.p  _  inf^  lk(u,J  _  <r£ (33) 
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where  S5j4  is  the  class  of  “statically  admissible”  stress-fields  given  by 


8":={»-  =  <t(u"h'+U'')  |  ii"®’  £  V**'  and  L(u""-)  =  0,  L(u',)  =  f}  (34) 

Note  that  (33)  defines  a  discrete  statically  admissible  patch-projection.  In  the 
computations  we  employed  the  same  sampling  points  as  the  ones  given  above  for 
the  ZZ-estimator.  It  is  also  possible  to  define  a  continuous  statically  admissible 
patch-projection  analogously  as  in  (28),  (29). 

h.  Least-squares  fit  for  the  displacements: 

Find  HgA  6  USA  such  that 

IK  -  =  JgL  IK  ~  (35) 

where  USA  is  the  class  of  “statically-admissible”  displacement-fields  given  by 

USA  :=  {u  =  uH”  +  up  |  u“H"  e  and  L(u“,^,)  =  0,  L(up)  =  f  j>  (36) 

The  corresponding  stress-field  &x  is  computed  from 

ax  :=  a( ufA)  (37) 

Here  {C/  denotes  the  set  of  sampling  points  for  the  displacement.  For  the 
quadrilateral  elements  we  used  the  mapped  (p+l)x(p-fl)  Lobat  to- Legendre 
points.  For  the  triangles  we  used  as  sampling  points  the  nodes  of  the  Lagrangian 
elements  of  degree  p.  It  is  also  possible  to  define  a  continuous  version  of  (35)  by 
employing  the  usual  L2-norm. 

After  a  recovered  stress  has  been  obtained  for  each  patch  from  the  above  recipes 
a  continuous  piecewise  polynomial  recovered  stress  is  obtained  over  the  entire  mesh 
by  combining  the  patch-projections  as  given  in  (30)-(32). 

Remark  4-2.  It  should  be  noted  that  the  cost  of  constructing  the  recovered-fields 
using  the  above  “statically-admissible”  recipes  is  always  less  than  the  corresponding 
cost  in  the  ZZ- recipes  because  we  minimize  over  a  smaller  set.  For  example,  in  the 
case  of  elasticity  and  cubic  elements  the  computation  of  arxz  for  a  patch  require 
the  inversion  of  an  27x27  matrix.  On  the  other  hand  the  calculation  of  the  &gA 
for  both  recipes  given  above  requires  only  the  inversion  of  a  12x12  matrix.  This 
difference  is  more  pronounced  for  higher-order  elements. 

Remark  4~3.  It  is  obvious  that  the  quality  of  a  recipe  can  be  improved  if  all  the 
known  properties  of  the  exact  solution  Me  employed.  Hence  it  is  not  surprising 
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that  the  quality  of  the  “statically-admissible”  recipes  is  better  than  the  quality  of 
the  ZZ-recipes  (see  Section  8.6,  below). 

Remark  4-4'  The  fact  that  it  is  possible  to  recover  fluxes  with  higher  accuracy,  than 
by  direct  computation,  makes  it  possible  to  construct  a-posteriori  error  estimates 
for  the  energy-norm  of  the  error. 

Remark  4-5.  We  do  not  address  the  higher-order  recovery  of  the  solution  (i.e.  the 
temperature  or  the  displacement- vector)  although  it  could  be  done  in  a  similar 
way. 


5  Outline  of  the  theoretical  setting 

The  theoretical  setting  will  be  outlined  below  for  the  class  of  locally  periodic 
meshes  which  are  defined  as  follows.  Let  0  <  H  <  H° ,  x°  =  (x°,  x£)  €  ft , 

5(x°,/f):={x  =  (x1,x2)  |  |xi-x?|<ir,  *  =  1,2}  (38) 

and  assume  H°  is  sufficiently  small  such  that  5(x°,  H°)  C  ft.  Further,  let  7  be  a 
set  of  multi-indices  (i,j) ,  x(‘j)  =  (xj‘j\x£j))  €  ft  and 

c  (xK  h)  :=  S(x<*\  h)  C  5(x°,  H) ,  (i,j)e  7  (39) 

be  the  set  of  the  A-cells  (or  cells)  which  cover  exactly  5(x°,  H)  as  is  for  example 
shown  in  Fig.  2a,  2b.  We  will  refer  to  5(x°,  H)  as  the  subdomain  of  periodicity  of 
the  mesh  centered  at  x°.  We  will  denote  by  c  :=  5(0,1)  the  unit-  (master-)  cell 
c,  the  h-cell  is  an  h-scaled  and  translated  master-cell. 

Let  T  be  a  mesh  of  triangles  or  squares  on  the  master-cell  (the  master-mesh) 
and  be  the  mesh  on  c(x^\  h)  which  is  the  scaled  and  translated  image  of 
T.  Let  Th  €  T  and  Th(x°,H)  be  the  restriction  of  Th  on  5(x°,  H)  and  the 
restriction  of  TA(x°,  on  c(x^\h).  We  assume  that  ( i,j )  e  7 

i.e.  Th(x°,H)  is  made  by  the  periodic  repetition  of  the  A-scaled  master  mesh.  The 
family  T  =  {Th }  is  a  class  of  meshes  which  are  locally-periodic  in  the  subdomain 
S(x°iH).  Outside  the  subdomain  5(x°,  H)  the  mesh  is  arbitrary;  it  could  have 
curved  elements,  refinements,  etc. 

Let  Q  be  a  polynomial  of  degree  (p  +  1)  defined  over  the  master-cell  c  and  let 
T  be  the  master-mesh.  Then  denote 

P  :=  Q  -  Gi”  (40) 

where  is  the  interpolant  of  degree  p  of  the  function  Q  defined  over  the  master- 
mesh  T  (for  which  h  =  1).  Any  polynomial  of  degree  p  on  an  element  rk  belongs  to 
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Si(rk)  and  hence  any  polynomial  of  degree  p  on  S(x°,  H)  belongs  to  S\(Th(xa,  H)). 
It  follows  that  p  defined  in  (40)  is  c-periodic  (this  can  be  shown  exactly  as  in  [54]) 
and 


p(l,x7)=  p(-l,x2),  \x2\  <  1 

(41a) 

p(i,,l)  =  p(ii,-l),  |x,|<l 

(416) 

Let 

FTper(c)  :=  ju  G  Hl(c)  |  u  satisfies  (41)  j 

(42) 

and 

^i,per(c)  :==  |u  e  #per(c)|  ulf  €  Si(f)  V  f  G  rj 

(43) 

Farther  let  zp  G  SpER(c)  such  that 

VrG  5?>per(c) 

(44a) 

and 

JAp~*p)  =  o 

(446) 

Note  that  the  function  zp  exists  and  is  uniquely  determined  (we  will  compute  it 
numerically  in  the  examples).  Let  us  also  define  ip  G  /f1(c)  by 


ip  :=  p—zp  =  Q  —  w  where  ib  :=  (jf*7  +  zp  (45) 

Let  ipK  G  HpER(c{x^\  h))  be  the  function  ip,  defined  above,  scaled  and  translated 
onto  the  cell  c{x^\  h)  of  the  mesh  in  S(x°,H)  i.e. 

(x),  ir(*)=  h,w(i) '  ,=1’2’  <46> 

where  x  =  i(x  —  x^),  x  G  c(x^\h).  It  is  easy  to  see  that  iph  can  be 
n 

periodically  extended  over  S(x°,Hl). 

In  [47]  we  proved  the  following  theorem  for  Poisson’s  equation  based  on  the 
theory  of  interior  estimates  (see  [53]): 

Theorem  1.  Let  Hx  <  H  <  H°  and  assume  that  the  following  assumptions  hold 
with 
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a  = 


6p+  1 


0=  p+  1  -  e, 


6p  ’ 

Assume  that  the  exact  solution  u  satisfies 


c  =  <70  = 


1 


6(6p  + 1) 


where  a  :=  (at, a2),  Dau  := 


<  K  <oo,  0  <  |ar|  <  p  +  2 


dx?  dx 2s  ’ 
i22  =  53  ao  >  0  where 


or |  :=  Qj  -f  a2,  and 
aa  :=  (X>“m)(x°) 


|o|=p+l 

Further  assume  that  the  mesh  Tk  is  such  that 


lleall£.»(s(*®,»1 »  ^  Ch0Ht ,  with  0>{p+  1)  -  c 
Moreover  assume  that  the  meshes  Th  in  S(x°,  H)  are  such  that 


(47) 


(48a) 


(486) 


(49) 


CXH°  <h<  C2H°  (50) 

Then  for  any  x  6  5(x°,  Hx ) 

!H=ltH+A^  (51) 

with  <r0  >  0,  |A|  <  1  and  C  independent  of  h. 

Theorem  2.  Let  the  assumptions  of  Theorem  1  hold;  then 

|  F(u)  -  f(«h)  |  =  |  F(Q)  -  F(w)  |  +  A CW+<  (52) 

with  Oq  >  0,  |A|  <  1  and  C  independent  of  h. 

Remark  5.1.  The  theorems  assume  that  the  mesh  is  periodic  and  that  the  solution  is 
smooth  in  a  small  subdomain  (i.e.  5(x°,  H))  in  the  interior  of  the  domain.  Outside 
the  subdomain,  we  assume  neither  periodicity  of  the  mesh  nor  smoothness  of  the 
solution.  The  solution  may  have  algebraic-type  singularities  at  a  finite  number  of 
comer  points  or  points  of  abrupt  change  in  the  type  of  boundary-condition.  Here  it 
is  only  assumed  that  the  pollution-error  in  a  shrinking  mesh-patch  (i.e.  TK(z°,  Hx )) 
in  the  interior  of  the  subdomain  is  controlled  (see  [51]).  This  implies  that  the  mesh 
has  been  adequately  refined  in  the  neighborhood  of  all  singular  points. 
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Remark  5.2.  If  we  further  assume  that 


*M:=  ||  |  *•(“-«■»)  |  ||t.(r)  >«’,  C  >  0  (53) 

theorem  2  implies  that:  A  point  xT  in  the  element  f  is  am  rj%- superconvergence 

point  for  the  averaging  T  if  and  only  if  7i(xT;  F\  T\  Q;  t,  T)  <  ——  where 

100 

H{x\  F;  T\  Q\  f,  f)  :=  njag  0(x;  F;f;Q;f,T)  (54) 

where 


S(*;JP;^;0;f,f) 


|  fWKx) i-  nam 

*(t) 


0 


with 

*  :=  map  |  F(V»)(*)  | 


if  #(f)  #  0 

(55) 

if  ^(f)  =  0 


(56) 


Here  Q  denote  the  class  of  all  (p  +  l)-degree  monomials  which  occur  in  Taylor 
series  expansions  of  functions  from  U  (see  also  Section  6  below). 

Remark  5.3.  Assumption  (53)  cam  be  realized  by  imposing  additional  restrictions 
on  the  vadues  aa  of  the  (p  +  l)-derivatives  of  the  solution  at  x0.  This  assumption 
is  reasonable  because  we  are  interested  in  a  sufficiently  large  class  of  solutions  U. 

Remark  5.4 •  The  proof  of  theorem  1  in  [47]  was  based  on  vau-ious  interior  estimates 
for  the  error  in  finite  element  approximations  of  Poisson’s  equation,  especiadly  the 
results  given  in  [53].  It  is  very  plausible  that  analogs  of  these  results  hold  for 
finite  element  approximations  of  the  elasticity  equations  and  more  general  elliptic- 
systems  because  the  maun  ideas  of  the  proofs  of  these  results  caurry  to  the  general 
case.  To  our  knowledge  the  precise  details  for  the  elasticity  equations  axe  not 
available  in  the  open  literature.  Nevertheless  we  will  assume  the  validity  of  the 
amalog  of  Theorem  1  for  the  equations  of  elasticity. 

Remark  5.5.  The  above  theoretical  setting  was  also  employed  to. assess  the  quadity 
of  a-posteriori  error  estimators  in  [54]  and  [55]. 
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6  The  method  of  freezing  the  periodicity 

Let  us  assume  that  for  a  given  locally  periodic  grid  with  corresponding  periodic 
master-mesh  T,  given  material  orthotropy  and  given  class  of  solutions  U  which  are 
smooth  in  the  interior  of  the  domain  (i.e.  the  solutions  which  satisfy  (48a))  we 
consider 


/  nd  P+1  x 

Q  ’=  W|  -  52<*kQk(*UX2),  Qk(xx,X2)  =  f  (57) 

the  class  of  (p+  l)-degree  monomials  which  occur  in  all  (p+  l)-degree  Taylor-series 
expansions  of  functions  from  U.  Here  Qk,  k  =  1,. . .  ,nd  denotes  a  set  of  linearly 
independent  monomials  which  form  a  basis  for  Q.  For  example,  let  us  assume 
that  U  is  the  class  of  solutions  UG  which  are  smooth  in  the  neighborhood  of  the 
subdomain  5(x°,  H);  in  this  case  we  may  choose 

Qki*i,  *a)  ==  x?2-kxk2-' ,  l<k<nd  =  p  +  2  (58) 

and  we  obtain  the  class  of  all  (p  +  l)-degree  monomials  QG. 

If  we  are  interested  only  in  the  sub-class  of  harmonic  solutions  UH  (i.e.  the 
solutions  in  U  which  satisfy  Au  =  0)  we  will  take  Q  as  the  two-dimensional  linear 
space  of  harmonic  monomials  of  degree  (pH-1)  denoted  by  QH,  namely, 

Qh  :=  \qh\Qh(x1,x2)  =  ^QfcQ"(x„x2)l  ,  (59) 

'  *=i  > 


Qi{xXix1)  =  'R.e(zp^x),  Q"(xi,x2)  =  Im(zp+1),  z  =  xx  +  ix2  .  (60) 

For  the  elasticity  problem  the  class  of  “harmonic”  (p+  l)-degree  monomials  Q“H" 
is  a  four-dimensional  linear  space.  Explicit  expressions  for  a  set  of  basis  monomial 
are  given  in  [50]. 

In  the  previous  Section  we  outlined  Theorem  2  which  states  that  we  can  obtain 
k  the  asymptotic  values  of  the  error  in  the  recovered  quantity  for  any  smooth  solution 
u  in  the  interior  of  a  periodic  mesh-subdomain  from  the  solution  of  the  periodic 
boundary-value  problem  (44),  obtained  using  the  master-mesh  T  over  the  master¬ 
cell  c,  with  data  obtained  from  the  local  (p  +  l)-degree  Taylor-series  expansion  Q 
of  the  exact  solution.  In  order  to  apply  the  results  of  the  theoretical  study  to  the 
practical  meshes,  for  which  the  mesh  is  not  locally  periodic  (like  for  example  the 
mesh  shown  in  Fig.  4a),  the  following  technique  of  freezing  the  periodicity  will  be 
employed: 
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1.  Let  Wq  be  an  interior  patch  of  elements  of  interest  (shown  shaded  gray  in 
Fig.  4a).  Define  the  s-layered  patch  of  elements  u>*,  s  >  1,  surrounding  the 
patch  Uq  by 

:=  U  ux  »  “x  '■=  U  T'  (61) 

X€N(t 0  X€N(  r>) 

where  iV(r')  denotes  the  set  of  the  vertices  of  element  r',  u;x  is  the  patch  of 
the  elements  connected  to  vertex  X.  The  patch  u>*  is  shown  with  thick  black 
perigram  in  Fig.  4a.  Note  that  the  patch  uft  with  its  surrounding  s-layers  is 
meant  as  a  pattern  which  may  be  repeated  in  other  places  in  the  grid. 

2.  Complete  the  patch  w*  to  a  periodic-grid  over  a  slightly  larger  square  periodic- 
cell  which  encloses  the  patch  as  shown  in  Fig.  4b.  The  periodic-cell  is  then 
scaled  and  translated  to  the  unit  master-cell  c. 

3.  Assume  that  the  mesh  in  the  neighborhood  of  element  r  is  made  from  the 
periodic  repetition  of  A-cells  obtained  from  the  master-cell  as  shown  in  Fig.  4c 
(by  h-scaling  and  translation)  and  let  h  tend  to  zero.  Based  on  the  results 
of  the  theoretical  study,  we  have 

fa?  *&,(W;*t,H)  =  K^Q^f.t)  (62) 

fa?  (63) 

where  R,,  A  denote  the  quantities  defined  in  the  master-mesh  in  the  master-cell. 

Remark  6.1.  The  above  limits  hold  for  the  locally  periodic  meshes  under  the  as¬ 
sumptions  of  the  theoretical  analysis  outlined  in  Section  5  (see  [47]  for  the  details). 
Hence  for  the  general  grids  the  limit  should  be  understood  for  the  mesh  which  is 
constructed  by  freezing  the  periodicity  (as  is,  for  example,  shown  in  Fig.  4c). 


7  Determination  of  the  77%-superconvergence 
quantities 

The  asymptotic  values  of  7f(x;  J-1',  Q;  f,  T)  for  a  recipe  T  for  a  class  of 

solutions  Q  can  be  determined  using  numerical  optimization.  Let  us  consider  a 
uniform  subdivision  of  the  element  r  into  subtriangles  with  vertices  at  the  set  of 
points  E  :=  {£*}*=!•  The  value  of  the  relative-error  function  for  the  class  Q  at 
the  point  £k  is 
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n-(L;  Fi  T\  Q\  f,  T)  :=  max  -=i - 3 - : - 100  (64) 

max  |E«,W,  )(£;)! 

j=l,...,np  ,_j 

where  }J=i  is  a  basis  of  Q  and  {iBi  denote  the  set  of  finite-element  solutions 
corresponding  to  the  Qi ’s  according  to  (45).  The  function  ^“(x;  F;  T\Q\r,  T)  can 
be  defined  for  any  point  x  €  f  by  using  linear  interpolation  in  the  subtriangles.  The 
7%-superconvergence  regions  for  the  recipe  T  in  the  element  t  can  be  approximated 
from 


^)(<3;^;f,f)«{x€f  |  H~(x;  F',T\  Q\r,T)  <  7/%|  (65) 

We  will  call  the  above  approach  the  direct  approach.  It  is  also  possible  to  use  a 
simplified  approach  which  avoids  the  use  of  numerical  optimization  at  every  point. 
We  define  the  function 


»-(*-,  F;T-,Q-T-,f)  :=~ 


(66) 


where 


Zs  :=  min 

“  a.- 


/  .max  |E«iW,)(£j)|^ 

j=l,...,np  i~i 


nd 


(E  a?) 


1=1 


(67) 


The  quantity  Zs  can  be  computed  using  numerical  optimization.  Let 
)(C;  T-,  T.  t)  :=  (x  6  r  |  «E(x;  F;  F;  <2;  f,  f )  <  ^  } 


(68) 


denote  the  approximate  regions  of  q%-superconvergence  for  the  class  of  solutions 
Q  obtained  by  the  simplified  approach. 

Remark  7.1.  Note  that 


R#„,«2;.F;t,T)  £  n^(Q-,F-,t,T) .  (69) 

Therefore  the  simplified  approach  results  to  a  conservative  estimate  for  the  rj%- 
superconvergence  regions. 
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The  function  VF^(x]U\  !FX\T^\  t,T),  defined  in  (24),  can  be  approximated  by 
computing  the  numerator  and  the  denominator  in  (24)  as  in  (59).  Then  approxi¬ 
mations  of  the  region  A(U;frl;F2\T,T),  defined  in  (25),  can  be  obtained. 

Remark  7.2.  The  functions  defined  in  (64)  and  (66)  depend  on  the  set  of  points  E. 
To  ensure  good  accuracy  in  the  approximation  of  the  r?%-superconvergence  regions 
a  sufficient  number  of  points  must  be  employed. 

Remark  7.S.  In  eq.  (64)  we  define  the  relative  error  by  taking  the  maximum  over 
the  class  of  solutions  Q  because  we  do  not  know  the  solution  a-priori.  Obviously  we 
have  to  use  all  the  information  about  the  solution  that  we  have  in  our  disposition. 
This  is  the  reason  that  in  (64)  we  compute  the  maximum  for  the  class  Q  (i.e.  the 
class  of  polynomials  of  degree  (p  +  1)  which  satisfy  the  differential  equation). 

8  A  model  numerical  study  of  recipes 

We  will  now  present  a  model  study  of  the  recipes  given  in  Section  4  based  on 
the  definition  of  superconvergence  given  in  Section  3.  We  will  compare  the  quality 
of  the  recipes  with  respect  to  direct  sampling  (i.e.  direct  evaluation  of  the  quality 
F(u)  from  the  finite  element  solution)  based  on  the  following  criteria: 

1.  The  gain  in  the  value  of  the  minimum  relative  error. 

For  a  given  recipe  we  will  define  the  minimum  value  of  the  maximum  relative 
error  for  the  class  of  solutions  U  in  the  element. 

Vmin(F;Jr;U;T,Th)  =  minn(x',F,Jr;U;T,Th)  (70) 

In  this  measure  we  will  say  that  recipe  is  superior  to  recipe  T7  in  the  element 
t  for  the  class  of  solutions  U,  if 

rfcm(F-,T'\U-,r,Tk)  <  (71) 

2.  Gain  in  the  T)  %-superconvergence  regions. 

For  a  given  ifecipe  T  and  given  q,  0  <  rj  <  100  we  define  the  percentage  area  of 
the  element  in  which  the  relative  error  in  the  recovered  quantity  is  less  than  rj%, 
namely 


T\  IA\  r,  Th ) 


100 


(72) 


Here  the  symbol  |  •  |  indicates  the  area  of  the  set.  Given  a  value  of  77,  we  will  say 
that  Tx  is  superior  to  T7  in  the  element  t,  for  the  class  of  solutions  U  if 


24 


r *<M(F-,r'-,U;r,Tk)  > 


(73) 


3.  Robustness  of  a  recipe  Tx  relative  to  the  recipe  T1  at  an  arbitrary  point  relative 
to  direct  sampling 

For  a  given  recipe  Tx  we  define  the  percentage  area  of  the  element  in  which  the 
recipe  is  less  robust  than  recipe  namely, 


100 


(74) 


In  this  measure  we  will  say  that  the  recipe  Tx  is  more  robust  than  the  recipe  T2 
in  the  element  r  if 


<  r^>t(U,r\r-,T,Th)  (75) 

We  will  also  define  the  maximum  value  of  the  relative  robustness  in  an  element, 
namely, 


T^^,T\r^r,Th)  =  maxVnu)W,T\T\T,Th)  (76) 

Below  we  will  investigate  the  quality  of  the  several  recipes  by  employing  the 
above  defined  measures. 

In  the  numerical  examples  we  will  employ  three  versions  of  the  ZZ-recipe,  and 
the  “staticaHy-admissible”  recipes  introduced  in  Section  4.  We  will  use  the  follow¬ 
ing  abbreviations: 

1.  ZZ-1 :  Discrete  ZZ  recipe,  Tzz~x ,  defined  in  (26)  with  c-kt  = 

2.  ZZ-2 :  Continuous  ZZ  recipe,  Tzz~% ,  defined  in  (28)  with  c-jW  = 

8.  ZZ-8:  Discrete  ZZ  recipe,  Tzz~z ,  defined  in  (26)  with  c =  6ik6jt. 

4.  SA-1 :  “Statically- admissible”  recipe,  TSA~X ,  based  on  least-squares  fit  for  the 
displacements  defined  in  (38)- (34). 

5.  SA-2:  “StaticaHy-admissible”  recipe,  FSA~2 ,  based  on  least-squares  fit  for  the 
displacements  defined  in  (35)- (36). 

6.  DIR:  Direct  sampling  recipe,  TD1R ,  which  means  that  the  solution  quantity 
is  evaluated  directly  from  the  finite-element  solution. 

We  will  employ  the  methodology  of  the  paper  in  the  study  of  these  recipes. 
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8.1  Influence  of  the  choice  of  the  solution  quantity  F(u) 
and  its  modulus  |  •  | 

In  this  Section  we  will  show  that  the  quality  of  a  recipe  depends  on  the  solution 
quantity  of  interest  and  the  modulus  employed.  We  considered  the  case  of  elasticity 
problem  in  which  one  could  be  interested  in  the  error  in  one  stress  component  or 
the  errors  in  all  the  stress  components.  We  employed  three  types  of  measures  and 
moduli  to  compute  the  errors  in  the  stress  components,  namely: 

(0  =  ^n(u)  and  |  <7n(u)  |  :=  |<ru(u)|. 

(ii)  F( u)  =  /i(<r0(u))  and  |/,(^(u))|  :=  |/1(<x0(u))|,  where  7,(o :=  <r„ 

+<^22  ^33* 

(iii)  ^(u)  =  and  |  ^(u)  |  =  J2(<ry(u)),  where  J2(crtj) 

where  cr'i}  :=  -  ^((T-,). 

We  computed  the  7/%-superconvergence  regions  in  the  elements  in  u)  (shown 
in  Fig.  4a)  for  the  stresses  computed  using  the  finite  element  solution  (i.e.  di¬ 
rect  sampling).  In  Figs.  6a,  6b,  6c  we  show  the  regions  of  r/%-superconvergence 

.JV)  for  F(a)  =  <r„(u),  /,(»;,(«))  and  respectively. 

The  values  of  TDIR\  e""";  r(  ,Tk)  and  Fd,r;  Q'"', t.  , Th)  for  the 

three  measures  used  in  Figs.  6a-6c  are  given  in  Table  2a.  It  can  be  seen  that  when 
F(u)  =  /i(o,j(u))  or  <7(J(u),  the  percentage  area  of  the  element  with  less  than  10% 
relative  error  is  small  or  zero.  We  also  computed  the  17%-superconvergence  regions 
for  the  stresses  recovered  using  the  recipe  Tzz~x .  In  Figs.  7a,  7b,  7c  we  give 
the  regions  ,Th)  for  F(u)  =  <rn(u),  /^(u))  and  ^(u), 

respectively.  In  Table  2b  we  give  the  values  of  ri*in(F]Pzz~1;Q‘‘fr-,Ti,Th)  and 
Tn<r,%(F-,Fzz~1-,  Q“fr";ri  ,Th)  for  the  three  measures  employed  in  Fig.  6a-6c.  It 
can  be  observed  that  when  the  three  different  measures  are  employed  different 
gains  in  the  value  of  the  minimum  relative  error  and  the  r/%-superconvergence  re¬ 
gions  are  obtained.  By  comparing  Tables  2a,  2b  we  observe  that  the  ZZ-recipe 
is  much  superior  them  the  direct  sampling  in  all  the  measures  i.e.  there  are  sig¬ 
nificant  gains  in  the  minimum  relative  error  and  the  r/%-superconvergence  regions 
especially  for  small  rj%.  For  example,  in  the  J2 -measure,  which  takes  into  account 
all  stress-components,  we  have  for  Fzz~l 


’liU**; G‘H’-,ri,Th)  <  2.5%,  >  85% 


while  for  direct  sampling 
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Similar  gains  are  also  observed  in  the  other  measures.  Note  however  that  the 
minimum  value  of  the  maximum  relative  error  in  each  element  obtained  from  the 
ZZ- recipes  and  the  direct  sampling  are  very  similar. 


8.2  Influence  of  the  choice  of  the  definition  of  the  ZZ- 
recipe 

We  now  give  an  example  to  show  how  the  methodology  can  be  employed 
to  compare  the  quality  of  different  recipes.  In  particular  we  considered  the  three 
versions  of  the  ZZ-recipe,  namely  fzz~x,  Fzz~2  and  Tzz~z . 

We  considered  the  elasticity  problem  and  the  elements  in  u>j}  (shown  in  Fig.  4b) 
assuming  biquadratic  elements  (p  =  2).  We  first  computed  the  r?%-superconvergence 
regions  for  stresses  which  are  recovered  using  JrZZ~1  and  Tzz~2.  In  Figs.  8a  and 
8b  we  show  the  regions  7?”* (u) ( ® ;  F;  t4 , TA)  for  T  —  Tzz~x  and  T  =  Tzz~2 , 

respectively.  In  Figs.  8c  and  8d  we  give  the  regions  F;  rf ,  Th)  for  T  = 

Tzz~x  and  T  =  T22~2 .  In  Table  3a  we  give  the  values  of  ^-n(cr-;  T\  Q“H  ]  r, ,  Th) 
and  QUfr]  ,  Th)  for  the  three  recipes  T22~x ,  Tzz~2  and  Jr2z~z.  It 

can  be  seen  that  the  quality  of  fzz~x  and  T2Z~Z,  using  this  measure,  is  almost  the 
same  whereas  the  quality  of  fzz~2  is  superior  to  the  other  two,  in  this  measure. 

We  also  computed  the  regions  in  the  elements  where  the  recipes  T22~x ,  fzz~2 
and  fzz~z  are  less  robust  than  the  direct  sampling.  In  Figs.  9a,  9b,  9c  we  show 
the  regions  •4^U(M)(CH" \F\FDIR',r,Th)  for  T  —  Tzz~x ,  JrZZ~2  and  Tzz~z,  re¬ 
spectively.  The  dark  regions  in  Figs.  9a-9c  are  the  regions  where  the  recipe  T 
is  less  robust  than  direct  sampling.  In  this  measure,  Tzz~z  appears  to  be  su¬ 
perior  to  J:ZZ~X  and  Tzz~2 .  However,  this  is  not  true  when  we  employ  F(u) 
=  <7-j  to  compute  the  regions  of  relative  robustness  A<Ti](QUH”;F;JrDIR-,T,Th) 
which  are  given  in  Figs.  9d,  9e,  9f  we  show  the  regions  for  T  —  T2Z~X,  Tzz~2 
and  ?2Z~Z.  In  Table  3b  we  report  the  values  of  TDIR]  r,  Th )  and 

Tv>i(ai ijQUirtJc;J:DIR;T,Th)  for  T  —  f22~x,  Fzz~2  and  Tzz~z,  respectively.  It 
can  be  seen  from  Figs.  9d,  9f  and  from  Table  3b  that  the  quality  of  J:ZZ~X  and 

fZZ-3 

is  about  the  same  while  the  Fzz~2  is  superior  to  the  other  two,  in  this 
measure. 


8.3  Influence  of  the  class  of  solutions  and  the  material- 
orthotropy 

In  this  Section  we  give  examples  to  show  that  the  quality  of  the  recipes  depends 
on  the  type  of  solution  (“harmonic”  or  general)  and  the  material  orthotropy.  We 
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considered  the  problem  of  orthotropic  heat  conduction  and  the  elements  in  u ))t 
(shown  in  Fig.  4b)  assuming  biquadratic  elements  (p  =  2).  We  computed  the 
Tf %- superconvergence  regions  'R%*u)(Q',Fzz~i;ri  ,Th)  and  the  regions  of  relative 
robustness  Aq.^{Q\!Fzz~x\TDIR\ r,  ,Th)  for  the  recipe  Tzz~x . 

We  first  considered  the  case  of  isotropic  heat  conduction  and  studied  the  qual¬ 
ity  of  the  recipe  for  the  class  of  harmonic  and  general  solutions.  In  Figs.  10a, 
10b  we  give  the  r?%-superconvergence  regions  Tzz~x\t{  ,Th)  for  Q  =  QH 


dxt 

and  Q  =  respectively.  Note  that  in  this  case  qAu)  =  — — ,  and 

VtT  | 


{( 


du 

dn 


du  ^  ^  du 


_  i  r  I  „  ,  .  It  can  be  seen  that  for  class  of  general  solutions  the 

dx\)  \dx2J 

superconvergence  regions  are  smaller  than  those  for  the  class  of  harmonic  solutions. 
In  Figs.  10c,  lOd  we  give  the  regions  of  relative  robustness  A ^  [Q,\TZZ~X  ,TDIR\ 

r,-  ,Th)  for  Q  =  QH  and  Q  =  QG,  respectively.  We  observe  that  the  regions  Ap 

{Q’,FZZ~X ,J:DIR\Ti  ,Th)  for  Q  =  QH  are  much  smaller  and  are  included  in  the 
corresponding  regions  for  Q  =  QG . 

We  now  consider  the  case  of  orthotropic  heat  conduction  with  the  principal  ma¬ 
terial  directions  oriented  along  the  coordinate  axes  and  ( — ^  )  =  2.  In  Figs.  11a, 

lib  we  show  the  regions  of  fj%-superconvergence  ( Q 5  JrZZ~ 1 ;  r,  ,Th)  for  Q  = 

QH  and  Q  =  QG,  respectively.  In  Figs.  11c,  lid  we  give  the  regions  of  relative  ro¬ 
bustness  -^(u)  ( £?;  JrZZ~x ,  TDIR\  r, ,  Th )  for  Q  =  QH  and  Q  =  QG,  respectively.  We 
reported  the  values  of  17*  n(<fc ;  Tzz~x ;  Q\  r- ;  Th)  and  T*<w%{q{ ;  Tzz~x\  Q;  t-jT,,)  in 
Table  4a  and  the  values  of  r^>1(9l ; Tzz~x  ^TD1R\  Q;  r,  -,Th)  in  Table  4b.  From  Ta¬ 
ble  4a,  it  can  be  seen  that  the  quality  of  the  JrZZ~x  recipe  deteriorates  when  the  ex¬ 
act  solution  is  not  “harmonic”  and  also  when  material  orthotropy  is  introduced.  In 
Table  4b  we  observe  that  regions  of  relative  robustness  *4?i(u)  ( Qi  Tzz~x ,  !FDIR\  r, ,  Th ) 
are  relatively  large  when  general  solutions  are  employed  or  when  the  material  is 
orthotropic.  Thus,  the  robustness  of  the  recipe  fzz~x  is  reduced  when  a  larger 
class  of  solutions  and  material-orthotropy  are  admitted. 

We  also  studied  the  influence  of  material  orthotropy  on  the  recipe  fzz~x  in 
meshes  of  triangular  elements  shown  in  Fig.  12.  Here  we  considered  the  peri¬ 
odic  mesh  of  triangles  shown  in  Fig.  12d  (Criss-Cross  Pattern)  with  p  —  2  and 
computed  the  values  of  »?mm(?«  >  CG;  r,  >  Th )  and  Tn<\o%  (?< ,  QG;  Ti ,  Th )  for 

T  =  J:ZZ~X  and  T  =  !FDIR.  These  are  given  in  Table  5  for  ^xx 


K 


K , 


=  1  and 


22 


-77—  =  10  with  the  principal  material  directions  oriented  along  the  coordinate- 
K  22 

axes.  In  Table  5  we  also  report  the  values  of  ( QG;  ?2Z~X ,  TD1R\  r,  ,Th)  and 
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T*>i(9ti  Qa\JrZZ~XiJrDIti\Ti ,  Th)  for  the  case  of  -  =  10.  We  observe  that  when 

K22 

the  ZZ  —  1  recipe  is  employed  there  is  a  significant  gain  in  the  minimum  relative 
error;  for  example  the  minimum  relative  error  obtained  with  direct  sampling  is 
greater  than  25%  while  that  of  the  ZZ  —  1  recipe  is  less  than  3%.  Moreover  there 
are  no  10%-superconvergence  regions  for  the  direct  recipe  while  there  are  regions 
of  significant  size  for  the  ZZ  —  1  recipe. 

8.4  Influence  of  the  polynomial  degree  of  the  elements 
and  the  topology  of  the  mesh  and  the  element-type 

We  now  give  some  examples  to  demonstrate  the  influence  of  element  type 
(quadrilateral  or  triangle,  the  degree  p  of  the  element)  and  the  local  mesh  topology. 

We  first  considered  the  elasticity  problem  and  the  elements  in  ujjj  (shown  in 
Fig.  4b).  In  Section  8.2,  we  used  this  mesh  patch  with  biquadratic  elements  to 
report  the  regions  of  relative  robustness  for  f22~x .  We  now  computed  these  regions 
for  p  =  1  and  p  —  3.  In  Figs.  13a,  13b  we  show  the  regions  Aa.} 
fDIR]Ti,Tk)  for  p  =  1  and  p  =  3,  respectively.  In  Table  6  we  give  the  values  of 
(Q*ir]r22-']jrDIn.TiiTh)  and  ,Th)  for  p  =  1, 

2  and  3.  It  can  be  observed  that  the  regions  of  relative  robustness 
jrzz~\ .  jrD/n.  T,  ^  7^)  where  the  ZZ— 1  recipe  is  less  robust  than  the  direct  sampling 
increase  in  size  as  p  is  increased. 

We  now  consider  the  problem  of  isotropic  heat  conduction  and  the  elements 
in  the  mesh  patch  shown  in  Fig.  4c  assuming  p  =  1.  We  first  computed  the 

q%-superconvergence  regions  of  the  derivative  — —  of  the  finite  element  solution 

OX  1 

(direct  sampling)  for  the  class  of  harmonic  and  genera]  solutions.  These  are  shown 
in  Figs.  14a  and  14b.  It  can  be  observed  that  the  q%-superconvergence  regions 
for  the  class  of  general  solutions  are  smaller  than  those  for  the  class  of  harmonic 

du 

solutions.  In  Figs.  15a,  15b  we  show  the  9  %- super  con  vergence  regions  of  — — 

OX\ 

which  is  recovered  using  the  recipe  J:ZZ~'1  for  the  class  of  harmonic  and  gen¬ 
eral  solutions,  respectively  We  also  computed  the  regions  of  relative  robustness 
Ag  ( Q\  Tzz~x\ TDlR\  r,- ,  Jh)  in  these  elements;  we  observe  that  the  regions  where 

the  ZZ—  1  recipe  is  less  robust  them  the  direct  sampling  increase  in  size  as  the  class 
of  solutions  is  enlarged.  In  Figs.  15c,  15d  we  give  these  regions  for  Q  =  QH  and  Q 
=  Qa,  respectively.  It  can  be  clearly  seen  that  r?%-superconvergence  regions  and 
the  regions  of  relative  robustness  depend  not  only  on  the  geometry  of  the  element 
but  also  on  the  topology  of  its  neighborhood. 

In  order  to  study  the  effect  of  mesh  topology  in  triangular  elements,  we  con¬ 
sidered  the  four  periodic  patterns  shown  in  Fig.  12  (Regular,  Chevron,  Union- Jack 
and  Criss-Cross  patterns).  We  solved  the  problem  of  isotropic  heat  conduction 
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using  these  patterns.  For  p  =  1,  the  ZZ-recipes  give  0%-superconvergent  deriva¬ 
tives  at  all  points  in  the  elements  of  Regular,  Union- Jack  and  Criss-Cross  pat¬ 
terns.  For  the  Chevron  pattern  and  p  =  1,  the  derivatives  recovered  using  the 
ZZ-recipe  are  0%-superconvergent  only  along  the  lines  which  pass  through  the 
midpoints  of  the  horizontal  sides.  In  Fig.  16  we  give  the  rj%- su percon vergence  re¬ 
gions  *£(<2°;*'ZZ-V,,:rfc)  in  the  elements  of  the  Chevron  pattern.  The  values 

_ m  ,  ,  %  ,  du 


of  1 1 


* 

min 


(fe;*,;QO;T4*Tjk)  “d  QG]Ti'Tk)  for  T  =  TZZ~'  “d  T 


—  7DIR  are  given  in  Table  7.  It  can  be  observed  that  for  7zz~l 

jrZZ'l;  T<  -T»>  >  43% 


whereas  for  direct  sampling 


%  <  &u 


Thus,  the  recovered  derivatives  in  the  elements  of  the  Chevron  pattern  are  not 
exact  at  all  points.  However  there  is  a  significant  gain  in  the  minimum  relative 
error  and  percentage  area  of  the  element  with  10%-relative  error. 

For  p  =  2,  the  recovered  derivatives  in  the  elements  in  all  the  four  patterns 
are  not  0%-superconvergent  except,  perhaps,  at  isolated  points  in  the  elements. 
In  Figs.  17a-d  we  give  the  >7%-superconvergence  regions  'RP^(QG,^zz~1\Ti , Th) 


for  the  four  patterns.  In  Table  8  we  give  the  values  of  >?^in(— — Qg;t,  ,Th) 

C/Xf 

and  r*<l0%(^-;7;  for  7  =  7ZZ~X  and  7  =  7DIR.  We  also  report 

\/Xf 


the  values  of  V$g)(&‘-,F2z-'-,FD'K-,Ti  ,T„)  and  ,Th) 

in  all  the  elements  of  the  four  periodic  patterns.  In  Fig.  17e  we  show  the  re- 
gions  'Rvg_(QG;  7zz~x,ri ,  Th)  for  the  distorted  Criss-Cross  pattern.  By  compar¬ 
ing  Figs.  17c  and  17e  we  observe  that  a  small  distortion  of  the  uniform  mesh  has 
caused  the  disappearance  of  5%-superconvergence  regions  in  one  of  the  elements. 


8.5  An  example  which  shows  that  the  results  obtained 
using  the  method  of  freezing  the  periodicity  hold  in 
the  actual  meshes 

We  now  give  some  examples  to  show  that  the  numerical  results  obtained  by  the 
method  of  freezing  the  periodicity  are  observed  in  actual  finite  element  calculations. 

We  first  give  am  example  to  show  that  the  pointwise  values  of  the  maximum 
relative  error  obtained  by  the  periodic-mesh  methodology  are  also  observed  in 
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actual  computations.  We  considered  the  problem  of  isotropic  heat  conduction 
and  studied  the  superconvergence  in  the  patch  ujj  assuming  biquadratic  elements 
( p  =  2).  We  used  the  methodology  of  Sections  6  and  7  to  compute  the  rj%- 
superconvergence  in  the  modulus  of  the  gradient  for  recipe  Tzz~x .  The  points  xfc 
where  the  maximum  relative  error  takes  a  minimum  value  in  the  element  are  shown 

in  Fig.  18,  i.e.  H(xk\  G“/r;  r{ ,Th)  =  Tzz~x\ f, ,  fh). 

The  values  of  the  maximum  relative  error  at  these  points  and  the  coefficients  of 
the  harmonic  monomials  which  correspond  to  the  extremal  values  of  the  error  at 
these  points  are  given  in  Table  9.  We  now  check  if  these  values  of  relative  errors 
can  be  observed  in  the  mesh  patch  when  solving  a  boundary-value  problem  in  the 
mesh  shown  in  Fig.  5a.  We  solved  a  Dirichlet  boundary  value  problem  with  the 
data  consistent  with  the  exact  solution 


ubx{x i,x2)  =  <>i(*i  ~  31,1a)  +  a2(I2  ~  3x?x2) 

where  the  values  of  a,  and  a2  are  given  in  Table  9.  The  values  of  a,  and  a3 
were  obtained  to  correspond  to  the  maximum  relative-error  at  the  point  xfc  in 
the  periodic  mesh.  The  values  of  the  relative  error  at  xk,  computed  using  the 
mesh  of  Fig.  5a,  are  given  in  Table  9.  It  can  be  observed  that  the  relative  errors 

at  these  points  are  close  to  the  values  of  H(xk;  \TZZ~X\ Q“ir\fi  ,Th)  (which 

are  computed  from  the  periodic-mesh).  Therefore  the  conclusions  made  from  the 
methodology  are  also  valid  for  finite  meshes. 

We  now  give  an  example  to  show  that  the  error  in  the  recovered  quantity 
rzz-\ 

(uk)  can  be  more  than  the  error  in  quantity  F(uh)  computed  directly  from 
the  finite  element  solution  at  some  points  and  for  some  solutions.  We  considered  the 
problem  of  isotropic  heat  conduction  and  studied  the  superconvergence  in  the  mesh¬ 
cell  Ujj  assuming  biquadratic  elements.  We  used  the  methodology  to  compute  the 
regions  fv  ,7\),  i.e.  the  regions  where  fz2~x  is  less  robust 

than  the  direct  sampling  TDlR.  We  chose  one  point  xk  in  each  element  in  u fjj 
which  lies  in  -4^.(2°; Tzz~x, TDlR\  ,Tfc)  (as  shown  in  Fig.  19)  and  determined 

the  coefficients  of  the  monomials  which  correspond  to  the  maximum  relative  error 

in  the  recovered  quantity  Fzz~x(uh)  at  this  point,  H(xk;  Fzz~x\  QG;  f,  ;T). 

(/X| 

These  coefficients  are  given  in  Table  10b.  We  then  solved  the  Dirichlet  boundary 
value  problems  using  the  mesh  shown  in  Fig.  5a.  The  data  for  each  problem  is 
consistent  with  the  exact  solution 

VEx(*l>*  2)  =  <*1*1  +  a2XlX2  +  a3XlX2  +  <*4*3 

where  the  values  of  {<*,  }£_,  are  given  in  Table  10b.  In  order  to  compute  the  error 
at  the  point  xk  we  chose  the  exact  solution  with  the  coefficients  of  the  monomials 
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corresponding  to  this  point  given  in  Table  10b.  The  values  of  the  errors  in  the 
solution  quantity  computed  directly  from  the  finite  element  solution  and  the  errors 
in  the  recovered  quantity  are  given  in  Table  10a.  We  also  computed  the  ratio 
of  the  error  in  the  recovered  quantity  j  ^zz~1(vh)  —  F(u)  |  to  the  error  in  the 
quantity  computed  directly  from  the  finite  element  solution  |  FDlR(uK)  —  F(u)  | ; 
these  values  are  given  in  the  last  column  of  Table  10a.  It  can  be  observed  that 
error  in  the  recovered  quantity  exceeds  the  error  in  the  finite  element  solution  by 
a  large  factor,  as  it  was  predicted  by  the  method  of  freezing  the  periodicity. 

8.6  The  quality  of  the  “statically  admissible”  recipes 

We  will  now  give  examples  to  illustrate  that  the  “statically-admissible”  recipes 
defined  in  Section  4.2  are  more  robust  than  the  ZZ-recipes  described  in  Section 
4.1. 

We  first  considered  Laplace’s  equation  and  the  elements  in  the  mesh-cell  u) 
(resp.  u>[[)  shown  in  Fig.  4a  (resp.  Fig.  4b)  assuming  bi-cubic  elements.  We 
computed  the  regions  (QH\  T\  rf- ,  Th)  for  the  levels  17  =  0.5%,  1%,  2%  for  T  = 

Tzz~ *,  PSA~l,  FSA~2.  These  regions  are  given  in  Fig.  20a,  20b,  20c  (resp.  Figs. 
21a,  21b,  21c).  We  note  that  for  the  recipe  Tzz~x  the  2%-superconvergence  regions 
are  very  small  or  non-existent  in  most  of  the  elements  while  for  the  recipes  PSA~l 
and  J:SA~2  the  2%-superconvergence  regions  cover  most  of  the  elements. 

We  also  considered  the  elasticity  problem  and  the  elements  in  the  mesh-cell 
ujjj  which  is  shown  in  Fig.  4c.  We  employed  biquadratic  elements  and  computed 
the  regions  «)  (G“ir;*;ii  ,Th)  with  t]%  =  0.1%,  0.5%,  1%  for  7ZZ~\  FSA~\ 
TSA~2\  these  regions  are  given  in  Fig.  22a,  b,  c.  We  observed  that  the  recipe  SA—2 
which  employs  the  least-squares  fit  of  the  displacement  is  the  most  robust.  For  this 
recipe  the  1%-superconvergence  regions  cover  almost  completely  all  the  elements 
in  Uf/j.  Note  also  that  both  “statically-admissible”  recipes  are  more  robust  than 
the  ZZ-recipe  at  all  the  points  in  all  the  elements. 

9  Summary  of  conclusions 

This  study  presented  a  new  definition  of  superconvergence  which  extends 
the  classical  definition  to  practical  meshes.  Using  this  definition  we  constructed  a 
methodology  for  checking  recipes  for  the  recovery  of  the  stresses,  fluxes  or  deriva¬ 
tives  from  the  finite  element  solutions.  We  employed  the  theory  to  assess  the  qual¬ 
ity  of  a  class  of  recipes  known  as  ZZ-recipes  and  improved  “statically-admissible” 
recipes.  The  following  conclusions  were  drawn: 

1.  The  regions  of  small  relative  error  (i.e.  less  than  10%)  for  the  ZZ-recipes  are 
much  larger  than  the  corresponding  regions  for  direct  sampling. 
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2.  The  performance  of  the  recipes  depends  on  the  class  of  meshes,  solutions  and 
materials  of  interest. 

3.  In  general,  we  cannot  say  that  a  particular  version  of  the  ZZ- recipe  is  the 
best. 

4.  Depending  on  the  geometry  of  the  mesh,  the  class  of  solutions  and  the 
material-orthotropy  significant  regions  of  large  relative  error  (greater  than 
50%)  may  exist  for  the  ZZ-recipes. 

5.  The  robustness  of  the  ZZ-recipe  relative  to  direct  sampling  decreases  when 
the  class  of  solution  is  enlarged  or  material-orthotropy  is  considered. 

6.  It  is  possible  to  improve  significantly  the  quality  of  the  recipes  by  employing 
the  available  information  about  the  classes  of  admissible  data.  Among  the 
recipes  tested  in  this  paper  the  recipe  SA-2,  which  employs  a  “harmonic” 
least-squares  fit  for  the  displacement,  is  the  best. 

7.  We  recommend  the  use  of  any  of  the  recipes  tested  here,  especially  the  recipe 
SA-2,  for  practical  computations. 

8.  The  conclusions  obtained  using  the  methodology  of  the  paper,  to  validate 
the  recipes,  hold  for  practical  computations. 

9.  The  methodology  should  be  used,  instead  of  random  benchmarks,  to  assess 
the  quality  any  concrete  recipe  used  in  practical  finite  element  codes  for  the 
classes  of  meshes,  solutions  and  materials  of  interest. 

10.  By  employing  the  methodology  of  freezing  the  periodicity  we  can  find  the 
regions  of  least-error  for  any  recipe  in  the  interior  of  any  practical  mesh  for 
the  various  solution  quantities  in  the  element. 

11.  The  presented  theory  has  asymptotic  character.  Hence  the  above  conclusions 
are  practically  valid  when  the  solution  is  smooth  and  the  mesh  is  sufficiently 
refined.  The  cases  of  unsmooth  solutions  and  coarse  meshes  will  be  addressed 
in  forthcoming  papers. 
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List  of  Figures 


Fig.  1.  Examples  of  geometrical  quantities  used  in  the  definition  of  the  recipes. 
The  value  of  the  recovered  quantity  at  the  point  x,  shown  here  in  the  interior  of 
mesh-cell  of  bilinear  elements,  is  obtained  as  a  linear  combination  of  the  values  of 
the  solution  at  the  vertices  {x}-®r 


Fig.  2.  An  example  to  explain  the  new  definition  of  superconvergence,  (a)  A 
typical  finite  element  mesh,  bilinear  elements,  in  the  sequence  of  meshes.  The 
point  (0.25,  0.75)  and  the  mesh-cell  consisting  of  four  sub-cells  (shown  shaded); 

(b)  The  mesh-cell  of  size  (2 h  x  2h)  and  the  x^' .  Here  we  show  the  points  Xj1 ,  x2; , 


Fig.  3.  An  example  to  explain  the  new  definition  of  superconvergence.  Exact  so¬ 
lution  uEX  (ij,  x2)  =  sin(xx1)sinh(jrx2);  bilinear  Elements.  Graphs  of  the  relative 
error  at  the  points  marked  in  Fig.  2b  for  various  values  of  the  mesh-cell  size  h.  Note 
that  the  values  of  the  relative  error  converge  to  the  theoretical  values  obtained  from 
the  periodic  mesh. 


Fig.  4.  The  interior  mesh-*  s  ti ployed  in  the  study  of  the  quality  of  the  recipes. 
(a)  Mesh-cell  ojj  consisting  of  mree  elements  shown  with  four  surrounding  mesh- 
layers;  (b)  Mesh-cell  u)jj  consisting  of  seven  elements  shown  with  seven  surrounding 
mesh-layers;  (c)  Mesh-cell  Ujjj  consisting  of  20  elements  surrounded  by  three  mesh- 
layers. 


Fig.  5.  The  method  of  freezing  the  periodicity,  (a)  Mesh  of  skewed  quadrilateral 
with  the  mesh-cell  Uq  shown  shaded  gray;  (b)  The  periodic  grid  obtained  by  com¬ 
pleting  the  mesh-cell  a/* ;  (c)  The  locally  periodic  mesh  obtained  by  the  method  of 
freezing  the  periodicity. 

Fig.  6.  Influence  of  the  choice  of  the  solution  quantity  F(u )  and  its  modulus  j  •  I 
on  the  quality  of  the  recipes:  Elasticity  problem,  biquadratic  elements,  “harmonic” 
solutions,  mesh-cell  u>j. 

(a)  The  regions  ***.,«•»•;  J***-,  r< ,  71); 

(b)  The  regions  r,  ,Tk),  /,( n)  :=  <rn(u)  +  <r22(u); 

(c)  The  regions  where  |cr„|  = 

The  levels  of  rj%  employed  are:  10%,  25%,  50%  (dark,  light  and  lighter  gray, 
respectively). 

Fig.  7.  Influence  of  the  choice  of  the  solution  quantity  F(u )  and  its  modulus  J  •  I 
on  the  quality  of  the  recipes:  Elasticity  problem,  biquadratic  elements,  “harmonic” 
solutions,  mesh-cell  Uj . 
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(a)  The  regions  IQ* (tt)(Q“H" ;  Tzz~x ;  t,  ,  TJ; 

(b)  The  regions  ft /*„)(£?“"" ; T22~' \ t,  , Th),  7,(u)  :=  <ru(u)  +  <r22(u); 

(c)  The  regions  K'l*{u)(QUfrlfzz-1\ r,  ,Th)  where  |  <rtJ  |  =  J2(cl}). 

The  levels  of  i)%  employed  are:  5%,  10%,  25%  (dark,  light,  lighter  gray,  respec¬ 
tively). 

Fig.  8.  Influence  of  the  choice  of  the  definition  of  the  ZZ-recipe :  r?%-superconver- 
gence  regions  for  the  ZZ-recipes.  Elasticity  problem,  biquadratic  elements,  “har¬ 
monic”  solutions,  mesh-cell  u>)j. 

(a)  The  regions  *£,„( 

(b)  The  regions  K*i.)(Q‘'r  ■,F2z-'\t>  ,Tk)-, 

(c)  The  regions  R;*ll)((2‘H';/’zz-1;Tj,Tj)  where  | <ry  |  =  djio,,), 

(d)  The  regions  where  |<r0  |  =  J,( <rtf). 

The  levels  of  rj%  employed  are:  5%,  10%,  25%  (dark,  light,  lighter  gray,  respec¬ 
tively). 

Fig.  9.  Influence  of  the  choice  of  the  definition  of  the  ZZ-recipe :  Regions  of  relative 
robustness  of  the  ZZ-recipes  with  respect  to  the  direct  recipe.  Elasticity  problem, 
biquadratic  elements,  “harmonic”  solutions,  mesh-cell  u>jj. 

(a)  The  regions  ,Tk); 

(b)  The  regions  ,Tk); 

(c)  The  regions  •A.u(.)(C*,r';  Fzz~3-,  FDm\  75  ,T»); 

(d)  The  regions  A.iM(Q'H'^z2-l\rmB-,Ti,Tk)-, 

(«)  The  regions  ;  ^zz  1; JFD,B; 7,  ,Tk); 

(f)  The  regions  A., w(a"'r;^zz-3;^°",;rj  ,Tt). 

The  dark-regions  indicate  the  areas  where  the  ZZ-recipes  are  less  robust  than  direct 
sampling.  In  the  regions  given  in  (d),  (e)  and  (f)  the  modulus  j  |  :=  J2(<rtj) 
was  employed. 

Fig.  10.  Influence  of  the  class  of  solutions :  Regions  of  77%-superconvergence  and 
of  relative  robustness  (with  respect  to  direct  sampling)  for  the  ZZ-recipes.  Isotropic 
heat-conduction,  biquadratic  elements,  mesh-cell  u ;*/. 

(a)  The  regions  R,i(C'';JrZZ-1;ri  ,Tk); 

(b)  The  regions  ,Tk)\ 

(c)  The  regions  A^_{QH\J:ZZ~^  ,TDtR\ T;  ,Tk); 

(d)  The  regions  A^{Qa^zz~\fDIR-,T,  ,7V). 

The  dark-regions  indicate  the  areas  where  the  ZZ-recipes  are  less  robust  than  direct 

du  \  2  /  du 

dx\ '  'cfcr2 


\  2 

J  was  employed. 


sampling.  Here  the  modulus 
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Fig.  11.  Influence  of  the  class  of  solutions :  Regions  of  r?%-superconvergence 
and  of  relative  robustness  (with  respect  to  direct  sampling)  for  the  ZZ-recipes. 

K 

Orthotropic  heat-conduction  with  =  2,  biquadratic  elements,  mesh-cell  u ))t. 

(a)  The  regions  *$.,( 

(b)  The  regions  RJ*  ,«2 ,Th)-, 

(c)  The  regions  AtM(QH\Fzz-',FDm\Ti  ,Tk); 

(d)  The  regions  A,iM(Qa;Fzz-\rD",-,Ti  ,Th). 

The  dark-regions  indicate  the  areas  where  the  ZZ-recipes  are  less  robust  than  direct 
sampling.  Here  the  modulus  |  q,  |  :=  \Jq\  +  q\  was  employed. 

Fig.  12.  Periodic  meshes  of  triangles,  (a)  Regular  pattern;  (b)  Union- Jack  pat¬ 
tern;  (c)  Chevron  pattern;  (d)  Criss-Cross  pattern. 

Fig.  13.  Influence  of  the  polynomial  degree  of  the  elements:  Elasticity  problem, 
tensor-product  quadrilaterals  “harmonic”;  solutions  mesh-cell  ujj.  The  regions 
A.l,(*f(Q',r\?zz-'\rD'R;Ti,Th). 

(a)  Bilinear  elements  (p  =  1); 

(b)  Bicubic  elements  (p  =  3). 

The  dark  regions  indicate  the  areas  where  the  ZZ  recipes  are  less  robust  than  direct 
sampling. 

Fig.  14.  Influence  of  the  topology  of  the  mesh  and  the  class  of  solutions:  Isotropic 
heat-conduction,  bilinear  elements,  direct  sampling. 

(a)  The  regions  R"Jl(<2";  ^“i  ’’i  >  T»)i 

(b)  The  regions  1V$L(Qai?D",;Ti  ,Th)\ 

The  levels  of  t)%  employed  are:  10%,  30%,  50%  (dark,  light,  lighter  gray,  respec¬ 
tively).  Note  the  17%-superconvergence  regions  for  the  general  class  of  solutions 
shown  in  (b)  are  smaller  than  the  corresponding  regions  for  the  harmonic  class  of 
solutions  shown  in  (a). 

Fig.  15.  Influence  of  the  topology  of  the  mesh  and  the  class  of  solutions:  Isotropic 
heat-conduction,  bilinear  elements,  mesh-cell  W///. 

(a)  The  regions  r, ,  T,); 

(b)  The  regions  ,Tt); 

The  levels  of  r\%  employed  are:  5%,  10%,  25%  (dark,  light,  lighter  gray,  respec¬ 
tively).  Note  in  contrast  to  the  ^%-superconvergence  regions  for  direct  sampling 
shown  in  Fig.  11  the  25%-region  practically  covers  all  the  elements. 

(c)  The  regions  ,Th)-, 
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(d)  The  regions  -dgL (QG;  Tzz~x ;  TDIR\  T{ ,  Th ). 

The  dark-regions  in  (c)  and  (d)  indicate  the  areas  where  the  ZZ-recipes  are  less 
robust  than  direct  sampling.  Note  that  the  dark-regions  for  the  general  class  of 
solutions  are  larger  than  the  corresponding  regions  for  the  harmonic  solutions. 

Fig.  16.  Influence  of  the  element  type :  Isotropic  heat-conduction,  linear  triangular 
elements,  Chevron-pattern,  discrete  ZZ-recipe. 


Fig.  17.  Influence  of  the  element  type:  Isotropic  heat-conduction,  quadratic  tri¬ 
angular  elements,  discrete  ZZ-recipe. 

(a)  The  regions  TV^{QG\Tzz~l\  r, ,  Th)  for  the  Chevron  pattern; 

(b)  The  regions  TV^{QG\ Fzz~l\Ti ,  Th)  for  the  Union-Jack  pattern; 

(c)  The  regions  TV^(QG\J:ZZ~l\ r, ,  Th)  for  the  Criss-Cross  pattern; 

(d)  The  regions  ( QG;  Tzz~x ;  t,  ,Ta)  for  the  Regular  pattern; 

(e)  The  regions  72.^ (QG;  Tzz~l\ri ,  Th)  for  the  distorted  Criss-Cross  pattern.  Note 

that  the  slight  displacement  of  the  central  node  caused  the  disappearance  of  the 
5%-regions  from  two  of  the  elements. 

The  levels  of  r]%  employed  are:  5%,  10%,  15%  (dark,  light,  lighter  gray,  respec¬ 
tively).  Here  the  modulus  j  ~~  j  :=  J (J“~)2  +  (J~)  was  empl°ye<^- 


Fig.  18.  An  example  to  show  that  the  conclusions  hold  for  general  meshes.  Mesh¬ 
cell  t*>jj  and  the  regions  TV^{Q.G\ Tzz~x ;  r{  ,Th)  for  biquadratic  elements  shown 

in  Fig.  10b.  The  points  of  interest  in  these  regions  are  marked  and  enumerated  for 
studying  the  relative  error  in  the  ZZ-recipe. 


Fig.  10.  An  example  to  show  that  the  conclusions  hold  for  general  meshes.  Mesh¬ 
cell  W//  and  the  regions  pzz~1\FDIR;Ti,Tk)  for  biquadratic  elements 

shown  in  Fig.  lOd.  The  points  of  interest  in  these  regions  are  marked  and  enumer¬ 
ated  for  studying  the  relative  error  in  the  ZZ-recipe. 


Fig.  20.  Comparison  of  the  robustness  of  the  “statically-admissible”  recipes  with 
the  ZZ-recipe:  Laplace’s  equation,  bi-cubic  elements,  mesh-cell  uj). 

(a)  The  regions  K'£.(QH;JrZZ-'-,T,,Thy, 

(b)  The  regions  ,Th); 

(c)  The  regions  n'Z(QH,FSA-,-,Ti  ,T,). 
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The  levels  of  r)%  employed  are:  0.5%,  1%,  2%  (dark,  light,  lighter  gray,  respec¬ 
tively).  Note  that  the  2%-superconvergence  regions  for  the  ZZ-1  recipe  are  almost 
non-existent. 

Fig.  21.  Comparison  of  the  robustness  of  the  “statically-admissible”  recipes  with 
the  ZZ-recipe :  Laplace’s  equation,  bi-cubic  elements,  mesh-cell  u ;*/. 

(a)  The  regions  VJ^{QH\JrZZ~x\Tx ,,TK)\ 

(b)  The  regions  (QH;  TSA~X ;  ); 

(c)  The  regions  7V^_  (QH]  FSA~2',  ). 

The  levels  of  r]%  employed  are:  0.5%,  1%,  2%  (dark,  light,  lighter  gray,  respec¬ 
tively).  Note  that  the  2%-superconvergence  regions  for  the  ZZ-1  recipe  do  not  exist 
in  four  elements  from  the  mesh-cell. 

Fig.  22.  Comparison  of  the  robustness  of  the  “statically-admissible”  recipes  with 
the  ZZ-recipe:  Laplace’s  equation,  bi-cubic  elements,  mesh-cell  ujjj. 

(a)  The  regions  ,Tk); 

(b)  The  regions  ;  r( ,  T»); 

(c)  The  regions 

The  levels  of  rj%  employed  are:  0.1%,  0.5%,  1%  (dark,  light,  lighter  gray,  respec¬ 
tively).  Note  that  the  1%-superconvergence  regions  for  the  SA-2  recipe  cover  the 
elements  almost  entirely. 
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Rate  of  convergence  of  the 
error  in  the  gradient 

Point  Xj  ;  9X  = 

19.7% 

h 

|f(«)-fZ2-'K)| 

h 

HH 

2.026 

8.104 

1.179 

9.432 

0.625 

10.001 

0.333 

10.640 

Point  x2  ;  02  =  20.92% 

h 

|F(u)-^2-‘K)| 

h 

0.25000 

1.917 

7.668 

0.12500 

1.143 

9.144 

0.613 

9.808 

HI 

0.317 

10.144 

Table  1.  Rate  of  convergence  of  the  error  in  the  gradient:  Laplace’s  equation, 
bilinear  elements  ( p  =  1),  meshes  consisting  of  mesh-cells  shown  in  Fig.2a.  The 
locations  of  the  points  xx, . . . ,  Xg  with  respect  to  the  mesh-cell  are  shown  in  Fig.  2b. 
Here  F(u)  denotes  the  value  of  the  derivative  of  the  exact  solution  Fzz~l(uh)  de¬ 
notes  the  recovered  derivative  using  the  recipe  described  in  Section  4  and  PDIR(uh) 
denotes  the  derivative  of  the  finite  dement  solution  uh. 
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Table  1.  (continued) 
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Table  1.  (continued) 
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7/%-superconvergence  for  the  ZZ- recipes 


Mesh-cell  of  seven  elements  shown  in  Fig.  4b 

Elasticity  problem 

■ 

Discrete- ZZ 

Continuous  ZZ 

Discrete- ZZ 

I 

with  aiikl 

(ZZ-  1) 

with  aiikt  {ZZ -2) 

with  6ik6t 

,  (ZZ  -  3) 

■ 

Vmin(ri  ) 

_% 

t-h<io% 

*l%in(Ti  ) 

T* 

Ttt<10* 

) 

TW<10% 

2.45 

11.96 

0.19 

61.37 

2.17 

12.94 

13.21 

■- 

2.74 

22.27 

13.62 

0.00 

14.10 

3.32 

37.16 

13.89 

0.00 

1.72 

1.61 

1.95 

42.22 

5 

2.16 

48.37 

3.11 

68.94 

2.09 

48.81 

□ 

9.31 

0.18 

7.59 

10.47 

9.33 

0.18 

H 

8.93 

3.27 

4.30 

39.31 

8.84 

3.33 

Table  3a.  Influence  of  the  definition  of  the  recipe:  Elasticity  problem,  “Harmonic” 
solutions,  biquadratic  elements  ( p  =2),  mesh-cell  Ujj.  rj  %-superconvergence  regions 
for  the  three  versions  of  the  ZZ  recipe.  The  values  of  ??*,„(  r, )  :=  T\  Q“H]  r, ,  Th) 

«<l  r*<1(M  ,T„)  for  T  =  Tz2-\  F**-*  and  Tz2-\ 
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Relative  robustness  of  the  ZZ-recipes 
with  respect  to  direct  sampling 

Mesh-cell  of  seven  elements  shown 

in  Fig.  4b 

Elasticity  problem 

■ 

Discrete- ZZ 

Continuous  ZZ 

Discrete- ZZ 

with  aukt  {ZZ  -  1) 

with  aiikt  {ZZ -2) 

with  Sik6it  {ZZ  —  3) 

■ 

TV>1 

ESI 

Tv>\ 

ES51 

Tv>\ 

i 

1.87 

4.81 

0.64 

0.00 

1.80 

4.79 

2 

2.31 

11.21 

1.79 

5.15 

2.09 

11.09 

3 

1.44 

9.44 

1.32 

3.20 

1.46 

9.52 

4 

1.39 

6.47 

0.83 

0.00 

1.36 

6.33 

5 

2.09 

10.39 

1.54 

5.61 

1.98 

10.34 

□ 

1.61 

3.01 

0.94 

0.00 

1.59 

3.01 

1 

2.42 

17.00 

1.05 

0.13 

2.42 

17.59 

Table  3b.  Influence  of  the  definition  of  the  recipe :  Elasticity  problem,  “Har¬ 
monic”  solutions,  biquadratic  elements  (p  =2),  mesh-cell  ujj.  Regions  of  relative 
robustness  of  the  ZZ-recipes  with  respect  to  direct  sampling.  The  values  of  Vpffi 


:=  max  DF(u)(x;  t{  ,Tk)  and  ,Th)  for  7  = 

7ZZ~\  7ZZ~2  and  7ZZ~3. 
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T7%-superconvergence  for  the  discrete  ZZ-recipe 


Mesh-cell  of  seven  elements  shown  in  Fig.  4b 


Orthotropic  heat-conduction 


1 

Q“ 

H” 

Q 

G 

H 

«21 

=  1 

9 

=  l 

Kn 

=  2 

■ 

K\i 

H 1 

Kn 

Kn 

1 

Vmin(Ti  ) 

T% 

TW<10% 

T* 

tH<10% 

Vm,n(Ti  ) 

T% 

Tw<io% 

Vmini'Ti  ) 

TW<10% 

1.86 

70.15 

n 

2.45 

11.96 

8.48 

0.09 

5.19 

24.47 

H9 

13.21 

0.00 

22.17 

0.00 

7.02 

14.01 

0.00 

14.10 

0.00 

19.56 

0.00 

1.72 

36.04 

1.66 

17.12 

1.72 

2.41 

pppa 

5 

2.22 

70.87 

3.48 

17.77 

2.16 

5.28 

mgm 

□ 

10.87 

0.00 

20.67 

0.00 

11.31 

0.00 

28.29 

0.00 

1 

4.78 

22.68 

9.70 

0.16 

8.93 

1.27 

11.97 

0.00 

Table  4a.  Influence  of  the  solution  type  and  material  orthotropy:  Orthotropic  heat 
conduction  problem,  biquadratic  elements  (p  =  2),  mesh-cell  a*;/.  »/%-superconvergence 
regions  for  the  discrete  ZZ  recipe  f22~x .  The  values  of  r/*in(rf )  :=  5  f2Z~x ;  Q;  r,-  ,  Th ) 

and  (q{ ;  Tzz~x\ Q\  rf ,  Tfc)  for  Q  =  Q“H”  and  Q  =  QG  with  the  principal 

K  K 

material  coefficients  =  1  and  -r—  =  2. 

K\\  Kn 


51 


Table  4b.'  Influence  of  the  solution-type  and  material  orthotropy :  Orthotropic 
heat  conduction  problem,  biquadratic  elements  (p  =  2),  mesh-cell  w)j.  Regions  of 
relative  robustness  of  the  ZZ-recipe  fzz~l .  The  values  of  t^>x  (qi ;  Q;  fzz~l }  J:DlR\  r, ,  Th ) 

for  Q  =  Q  and  Q  =  QG  with  material  orthotropy  =  1  and  =  2. 

Ku  I'll 
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»7%-8uperconvergence  and  relative  robustness  of  the  ZZ-recipe 
with  respect  to  direct  sampling 


Quadratic  triangles  (p=2) 

Criss-Cross  pattern 

j:ZZ-l 

■pDIR 

Kn 

*22 

*u  =  ,„ 

Kii 

*n  _  j 
*22 

|u  =  1„ 

A  22 

0 

A  22 

% 

Tn<  io% 

Vmin(Ti) 

r% 

TH<10% 

Vmin(Ti  ) 

T% 

tH<  10% 

T% 

'n<  io% 

x  100 

T% 

rZ»i 

2.86 

2.99 

2.86 

2.99 

2.58 

2.65 

2.58 

2.65 

1.96 

8.94 

1.96 

8.94 

1 

0.00 

0.00 

0.00 

0.00 

25.87 

24.55 

25.87 

24.55 

0.00 

0.00 

0.00 

0.00 

144.90 

173.78 

144.90 

173.78 

6.08 

7.59 

6.08 

7.59 

Table  5.  Influence  of  material-orthotropy:  Criss-Cross  pattern.  Orthotropic  heat- 
condiction,  quadratic  triangles  (p  =  2).  Values  of  T?*tn(T« )  a°d  Tn<io%  f°r  the  class 
of  general  solutions  and  the  recipe  fzz~l  (columns  1-4)  and  the  recipe  TDIR 
(columns  5-8).  In  columns  9-10  we  give  the  values  of  and  for  the  case 
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Relative  robustness  of  fZ2~x 
with  respect  to  TD1R 


Mesh-cell  of  seven  elements 
shown  in  Fig.  4b 

Elasticity  problem 

warn 

■a 

msm 

mm 

mm 

mm 

mm 

mm 

I 

1 

■ 

0.00 

0.00 

0.00 

0.00 

0.00 

0.00 

0.00 

1.14 

2.28 

2.12 

1.50 

2.95 

1.54 

1.36 

9.37 

15.53 
19.04 
13.64 
6.24 
20.98 

10.53 

3.34 

2.92 

5.65 

3.64 

3.38 

4.82 

2.90 

42.36 

79.31 

82.91 
33.17 
20.97 

67.92 
29.76 

Table  6.  Influence  of  the  p-order  of  the  element:  Elasticity  problem,  bi-p  elements, 
mesh-cell  u>jj.  Regions  of  relative  robustness  of  the  ZZ-recipe  J^zz~l  with  respect  to 
direct  sampling.  The  values  of  Vpffi  :=  max  X?F(uj(x;  Q “rr ;  J:ZZ~X , TDIR\  r, ,  Th) 

and  t*>x  {q{ ;  Q;JrZZ~1,JrDIR;  r{ ,Th)  for  p  =  1,  2  and  3. 
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t/%-superconvergence  and  relative  robustness 
of  the  ZZ-recipe  with  respect  to  direct  sampling. 


Isotropic  heat-conduction 

Linear  triangles  (p  =  1) 

Chevron  pattern 

JpZZ-l 

jtDIR 

Ti 

V minimi  ) 

T»<10% 

Vmin(Ti  ) 

T% 

~H<  10% 

r% 

TV>1 

1 

0.00 

43.97 

43.68 

0.00 

23.20 

0.00 

2 

0.00 

43.97 

43.68 

0.00 

23.20 

0.00 

3 

43.97 

41.91 

0.00 

22.59 

0.00 

4 

B 

43.97 

41.91 

0.00 

22.59 

0.00 

Table  7.  Influence  of  the  topology  of  the  mesh:  Chevron  pattern.  Isotropic 
heat-conduction,  linear  triangles  (p  =  1).  Values  of  the  minimum  relative  er- 

ror  Vmin(Ti )  :==  (  o — <2G;  Ti  i  Th)  in  the  elements  for  the  class  of  general  so- 

OXi 

dtx 

lutions  and  percentage  area  of  the  element  T«<io%  2G;t,  ,Ta)  in  which 

the  relative  error  in  the  recovered  gradient  is  less  than  10%  for  T  —  fzz~l 
and  T  =  TDIR  (columns  3,  4).  In  columns  5,  6  we  give  the  values  := 

max  (x;  QG\Tzz~l ,Tdir-,  t,  ,  Th )  and  of  t$>x  defined  in  (74). 
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T/%-superconvergence  and  relative  robustness 
of  the  ZZ-recipe  with  respect  to  direct  sampling. 


Isotropic  heat  conduction 


Quadratic  triangles  ( p  =  2) 


jrZZ-  1 


jrDIR 


Vrnin(Ti )  r7i<l0%  l?min(Ti )  TW<10%  xliwl  7"z?>1 


Regular  Pattern 


100.00 

100.00 


Chevron  Pattern 


0. 


29.34 

30.09 

61.72 

29.10 

61.72 

29.10 

Union  Jack  Pattern 


.00 


31.77 

31.77 

28.51 

28.51 

28.51 

28.51 

31.77 

31.77 


Criss-Cross  Pattern 


.35 
.35 
41.49 
41.49 


[iTlIll 


112.83 

112.83 

131.51 

131.51 


31.24 


134.29 


Table  8.  Influence  of  the  mesh-topology:  Regular,  Chevron,  Union-Jack  and  Criss- 
Cross  patterns.  Isotropic  heat-conduction,  quadratic  triangles  (p  =  2).  Values  of 
Vmin(Ti )  and  Tn<io%  f°r  ^e  re°jpes  T2Z~X  (columns  1,  2)  and  TDIn  (columns  3, 
4).  In  columns  5,  6  we  give  the  values  of  and  t*>x. 
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Coefficients  of  the  cubic  harmonic 
solutions  used  in  Table  9a 


u(a?,,x2)  =  Oj(xi  —  Zx^l)  +  a2(:r2  —  Zxjx3) 


Point  k 

mBm 

Or2 

1 

1.000 

0.714 

2 

1.000 

0.688 

3 

-0.577 

1.000 

4 

-1.000 

-0.796 

5 

-1.000 

-0.796 

6 

-0.578 

.1.000 

7 

-1.000 

-0.578 

Table  9b.  An  example  which  shows  that  the  results  obtained  using  the  method 
of  freezing  the  periodicity  hold  in  the  actual  meshes.  Isotropic  heat-conduction, 
biquadratic  elements  ( p  as  2),  mesh-cell  Ujj.  The  coefficients  of  the  solutions 
employed  to  compute  the  errors  at  the  points  xt  in  table  9a  (by  solving  a  Dirichlet 
boundary- value  problem  over  the  square-domain  ft  =  (0,  l)2  using  the  mesh  shown 
in  Fig.  5a). 
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Coefficients  for  the  cubic  solutions  used  in  Table  10a 


u(xl5x2)  = 

°lIl  +  <*2*1*2 

+  a3XjX3  -|-  a4 

X3 

*2 

Point 

<*i 

<*2 

<*3 

<*4 

1 

1 .00000000 

.27931142 

-.88763169 

.92384722 

2 

-.05439094 

-.02646028 

.84947474 

-1.00000000 

3 

.43162747 

-1.00000000 

-.93122513 

.45706703 

4 

-.75131974 

-.61557397 

.42834683 

1.00000000 

5 

.95275077 

-.75780325 

-1.00000000 

-.79558135 

6 

7 

-.24876338 

-.29648689 

.27911873 

-1.00000000 

1.00000000 

-.21416015 

.59140203 

.97156428 

Thble  10b.  An  example  which  shows  that  the  results  obtained  using  the  method 
of  freezing  the  periodicity  hold  in  the  actual  meshes :  Isotropic  heat  conduction, 
biquadratic  elements  (p  =  2),  mesh-cell  The  coefficients  of  the  solutions 
employed  to  compute  the  errors  at  the  points  xk  in  Table  10a  (by  solving  a  Dirichlet 
boundary  value  problem  over  the  square-domain  ft  =  (0,  l)3  using  the  mesh  shown 
in  Fig.  5a).  For  example,  to  compute  the  errors  at  point  xi  in  Table  10a  we 
employed  the  exact  solution  u(x1,x2)  =  x?  +  0. 27931 142x3xa  -  0.88763169x,x2  + 
0.92384722x2. 
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The  Laboratory  for  Numerical  Analysis  is  an  integral  part  of  the  Institute  for  Physical 
Science  and  Technology  of  the  University  of  Maryland,  under  the  general  administration  of  the 
Director,  Institute  for  Physical  Science  and  Technology.  It  has  the  following  goals: 

To  conduct  research  in  the  mathematical  theory  and  computational  implementation  of 
numerical  analysis  and  related  topics,  with  emphasis  on  the  numerical  treatment  of 
linear  and  nonlinear  differential  equations  and  problems  in  linear  and  nonlinear  algebra. 

To  help  bridge  gaps  between  computational  directions  in  engineering,  physics,  etc.,  and 
those  in  the  mathematical  community. 

To  provide  a  limited  consulting  service  in  all  areas  of  numerical  mathematics  to  the 
University  as  a  whole,  and  also  to  government  agencies  and  industries  in  the  State  of 
Maryland  and  the  Washington  Metropolitan  area. 

To  assist  with  the  education  of  numerical  analysts,  especially  at  the  postdoctoral  level, 
in  conjunction  with  the  Interdisciplinary  Applied  Mathematics  Program  and  the 
programs  of  the  Mathematics  and  Computer  Science  Departments.  This  includes  active 
collaboration  with  government  agencies  such  as  the  National  Institute  of  Standards  and 
Technology. 

To  be  an  international  center  of  study  and  research  for  foreign  students  in  numerical 
mathematics  who  are  supported  by  foreign  governments  or  exchange  agencies 
(Fulbright,  etc.). 

Further  information  may  be  obtained  from  Professor  I.  BabuSka,  Chairman,  Laboratory  for 
Numerical  Analysis,  Institute  for  Physical  Science  and  Technology,  University  of  Maryland,  College 
Park,  Maryland  20742-2431. 


